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PROJECT  SUMMARY 


This  report  studies  some  aspects  of  analysis  and  controller  design  for  linear  systems 
with  parameter  uncertainties.  It  is  assumed  that  the  uncertainties  are  quantified  as 
intervals  and  the  relevant  parameter  can  assume  any  value  from  the  interval. 

Unlike  conventional  methods,  use  of  interval  arithmetic  is  investigated.  Interval 
arithmetic  provides  a  mechanism  for  dealing  with  the  uncertainties  in  a  very  general 
framework.  All  the  calculations  have  to  be  performed  in  terms  of  intervals  instead  of 
real  and/or  complex  numbers. 

Historically,  the  use  of  interval  arithmetic  was  motivated  by  an  elfort  to  find 
bounds  on  round-off  errors  in  numerical  computations,  and  therefore,  typically  the 
interrals  used  were  of  the  order  of  the  last  significant  place  after  the  decimal  on 
the  computer  on  which  it  was  being  employed.  The  premises  in  this  report  is  that 
by  increasing  the  interval  size  to  accommodate  model  uncertainty,  if  a  controller  is 
designed  for  an  interval  instead  of  a  particular  value  from  that  interval,  stabili(y  of 
the  system  in  face  of  complete  variation  in  corresponding  parameter  is  autr  niatically 
guaranteed! 

At  the  current  time  the  above  objective  is  achieved  only  at  the  cost  of  overbound¬ 
ing  the  uncertainties  hence  the  controller,  when  it  exists,  will  be  conservative. 


IV 


1.  INTRODUCTION  AND  OVERVIEW 

1.1.  Introduction 

The  fields  of  optimization  and  control  theory  are  well  developed  disciplines  with  very 
strong  theoretical  and  methodological  foundations.  Although  a  significant  effort  has 
been  underway  to  develop  both  of  the  above  areas  independently,  a  large  variety  of 
current  problems  require  us  to  pool  the  resources  from  these  and  other  related  areas 
and  to  put  forth  a  combined  effort  to  solve  them.  Among  the  problems  presenting 
these  challenges  is  the  development  of  a  symbiotic  relationship  between  system  theory 
principles  and  practical  problems  such  as  design  of  flight  control  systems,  communi¬ 
cation  systems,  transportation  systems,  etc. 

Control  theory  has  been  a  subject  of  intensive  research  for  several  decades.  How¬ 
ever,  until  recently,  the  principle  effort  of  research  was  focussed  on  classical  design 
methodologies  such  as  Bode  diagram,  root  locus,  Nyquist  plots,  etc.  These  design 
methods  are  extremely  effective  for  solving  small  scale  and  relatively  less  complicated 
problems.  However,  they  have  limited  utility  in  modern  multivariable  control  sys¬ 
tems  analysis  and  design  due  to  the  stringent  design  and  performance  requirements 
commanded  by  current  technology. 

An  important  aspect  of  analysis  and  design  of  a  system  is  to  make  it  “robust” 
in  face  of  uncertainties  in  the  model  parameters  and/or  operating  conditions  [3],  [5], 
[7],  [9],  [12],  [16],  [35],  [40]  and  [44],  The  need  to  incorporate  robustness  in  design 
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is  necessitated  by  the  fact  that  for  most  practical  systems,  the  model  is  known  oiily 
approximately.  For  example,  in  the  aircraft  industry,  the  aircraft  model  i.s  constiucteci 
using  the  data  obtained  from  the  wind-tunnel  experiments  on  the  aircraft  body.  As 
a  consequence,  the  parameters  of  the  model  would  not  have  a  specific  value,  ratln'i 
they  are  known  to  lie  within  an  interval.  Since  the  actual  flight  data  are  not  availabk'. 
the  controller  should  be  able  to  account  for  the  unmodeled  parameters  that  can  l)e 
obtained  only  when  the  aircraft  is  airborne.  Design  of  a  ‘“robust"’  controller  therefori' 
becomes  a  priority  in  such  applications. 

Fortunately,  in  the  last  decade  the  momentum  has  visibly  shifted  to  integrating 
more  sophisticated  mathematical  theories  to  solve  these  problems.  As  a  consequence, 
the  newer  design  techniques  are  capable  of  meeting  the  required  performance  measures 
and  at  the  same  time  giving  the  designer  sufficient  freedom  to  incorporate  additional 
“desirable  features”  of  robustness.  Several  theories,  notably  those  developed  in  [3]. 
[9],  [11],  [13],  [14],  [19],  [24],  [33]  and  [44]  have  been  put  forth  that  enable  us  to  design 
robust  systems.  One  assumption  underlying  all  of  these  theories  is  that  the  normnal 
model  of  the  system  is  known.  The  design  is  carried  out  based  on  this  model,  and  the 
final  product  is  analyzed  for  the  degree  of  robustness. 

The  approach  adopted  in  the  propos«,'d  research  is  significantly  different  from  th( 
conventional  approaches  to  the  solution  of  this  problem.  In  this  project,  we  will  inves 
tigate  the  application  of  interval  mathematics  in  analysis  and  design  of  multivarial)Ie 
control  systems.  Interval  analysis,  as  introduced  by  Moon-  [26],  [27]  and  ('xtinided 
by  several  researchers  e.g.,  [1],  [30],  considers  each  nuinlier  as  an  int('rval  instt-ad  of 
a  fixed  point  in  the  conii)lex  plane.  .Although  the  original  motivation  Ix'hind  th<'  use- 
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of  interval  analysis  was  to  capture  all  rounding-errors  during  computations,  the  in¬ 
terval  arithmetic  was  developed  as  a  parallel  to  conventional  arithmetic.  Since,  as 
discussed,  the  parameters  of  the  models  in  modern  multivariable  systems  are  usually 
known  to  be  within  an  interval,  the  use  of  interval  arithmetic  for  their  analysis  and 
design  appears  to  be  a  logical  direction. 

The  prediction  and  control  of  dynamical  processes  described  by  general  physical 
laws  or  experimental  data  is  a  basic  problem  in  engineering.  Therefore,  the  starting 
point  in  analysis  and  design  of  a  practical  system  is  to  obtain  an  accurate  description 
or  a  model.  These  models  are  usually  obtained  by  means  of  simulations  and/or 
laboratory  experiments.  If  the  model  parameters  are  known  precisely,  the  solution  of 
this  representation  problem  is  very  well  understood.  In  practice,  the  real  conditions 
under  which  the  system  would  operate  can  only  be  approximated  in  the  laboratory 
experiments.  Hence  it  is  reasonable  to  assume  that  a  precise  knowledge  of  the  process 
or  model  parameters  is  atypical.  This  is  especially  true  in  the  case  of  the  aircraft 
industry  where  wind-tunnel  tests  are  used  to  obtain  the  best  possible  approximation 
of  the  “true”  model  of  the  aircraft.  The  control  elements  such  as  rudder,  elevator, 
etc.,  are  instedled  based  on  the  available  model. 

For  applications  such  as  the  one  described  above,  although  a  precise  model  is 
almost  never  available,  the  extensive  experience  with  such  experiments  over  the  past 
several  decades  does  enable  the  engineers  to  specify  the  degree  of  uncertainty  associ¬ 
ated  with  various  parameters  that  determine  the  performance  of  given  system. 

In  classical  frequency  domain  techniques  for  single-input,  single-output  systems, 
the  concept  of  gain  and  phase  margins  are  well  understood.  These  margins  provide 
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effective  measure  in  determining  the  robustness  of  the  given  plant.  These  notions, 
however,  do  not  directly  extend  to  multi-input,  multi-output  plants.  In  fact,  it  h:is 
been  shown  for  the  multivariable  case  that  even  if  each  channel  in  the  feedback  mech¬ 
anism  possesses  some  desired  gain  and  phase  margins,  a  small  perturbation  in  the 
system  can  result  in  closed-loop  instability. 

Analysis  and  design  of  a  system  whose  parameters  are  known  to  lie  within  a 
range  (rather  than  having  an  exact  value)  has  received  considerable  attention  over 
the  last  decade.  Various  analysis  and  design  methodologies  used  by  control  systems 
engineers  are  essentially  meant  for  application  to  a  “nominal”  model.  The  resulting 
design  is  said  to  be  robust  if  the  system  performs  within  acceptable  limits  in  the  face 
of  significant  parameters  variations  and  model  uncertainties. 

This  situation  has  led  to  fundamental  extensions  and  re-evaluations  of  the  mod¬ 
eling  philosophy,  an  extensive  and  rigorous  development  of  powerful  estimation  tech¬ 
niques  and  the  advent  of  various  approaches  to  the  design  of  controllers  that  are 
robust  ageunst  such  uncertainties.  Some  of  the  more  important  of  the  advances  that 
address  the  design  of  robust  controllers  are  (1)  Robust  Servomechanism  Control,  where 
a  good  understanding  of  the  nature  of  uncertainties  associated  with  the  system  and 
the  environment  in  which  it  operates  enables  us  to  design  servomechanism  controllers 
that  guarantee  asymptotic  regulation,  (2)  Linear  Quadratic.  Regulator  with  multiple' 
loop  state  feedback.  These  regulators  have  excellent  robustness  properties  when  inea 
sured  by  classical  criteria  of  gain  and  ])hase  margins  and  can  undergo  substantial  gain 
and  phase  perturbations  without  becoming  ur<^table,  and  (3)  Optimal  Hankel  Norm 
Optimization  based  techniques  that  use  matrix  valued  interpolation  theory  to  ol)tain 
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robust  controllers  for  the  given  plant. 

In  practically  all  of  the  research  that  has  been  conducted  in  robust  controller 
design,  the  underlying  model  is  assumed  to  have  some  nominal  values  of  the  param¬ 
eters.  This  choice  of  parameters  is  carried  through  the  design.  Once  a  controller  is 
found,  the  robustness  of  the  resulting  controller  is  determined  from  the  final  design. 
Intuitively,  it  appears  more  appealing  to  perform  the  complete  design  using  intervals 
in  which  these  parameters  lie,  instead  of  selecting  one  particular  value  from  these  in¬ 
tervals.  This  leads  to  the  natural  choice  of  interval  arithmetic  as  the  most  logical  tool 
to  address  the  problem  of  robust  controller  design.  The  subject  of  interval  analysis 
was  developed  in  mid-sixties  in  quest  for  rigor  in  numerical  computations  on  com¬ 
puters  using  finite  precision  arithmetic.  Interval  arithmetic  treats  intervals  as  a  new 
kind  of  number.  Computations  in  appropriately  rounded  arithmetic  produce  results 
that  contain  both  ordinary  machine  arithmetic  results  as  well  as  infinite  precision 
arithmetic  results. 

The  use  of  interval  arithmetic  was  motivated  by  an  effort  to  find  bounds  on  round¬ 
off  errors  in  numerical  computations,  and  therefore,  typically  the  intervals  used  were 
of  the  order  of  the  last  significant  place  after  the  decimal  on  the  computer  on  which  it 
was  being  employed.  However,  since  the  underlying  principles  for  various  operations 
such  as  addition  and  multiplication  of  the  intervals  have  already  been  established, 
it  remains  to  increase  the  interval  size  that  suits  the  practical  problem  of  enclosing 
the  uncertainty  in  the  model  parameters,  and  extend  the  methods  of  analysis  and 
design  of  systems  using  interval  arithmetic.  It  is  not  too  difficult  to  visualize  that  if  a 
controller  is  designed  for  an  interval  instead  of  a  particular  value  from  that  interval, 
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stability  of  the  system  in  face  of  complete  variation  in  corresponding  parameter  is 
automatically  guaranteed!  It  should  be  emphasized  that  this  is  achieved  only  at  the 
cost  of  overbounding  the  uncertainties  hence  the  controller,  when  it  exists,  will  be 
conservative. 

In  the  rest  of  this  report,  we  will  consider  system  described  by  their  state-sjjace 
representation: 

=  Ax(t)  +  Bu(t)  (1.1a) 

at 

y(t)  =  Cx(t),  (1.16) 

where  x(t)  G  IR",  u(t)  G  IR"^  and  y(t)  G  IR^.  Or  in  their  transfer  function  represen¬ 
tation: 


Y(s)  =  G(s)V{s) 

(1.2a) 

= 

dis) 

(1.26) 

=  C(sl  -  A)"^BU(s) 

(1.2c) 

In  most  existing  analysis  and  design  techniques,  it  is  assumed  that  the  matrices 
A,  B  and  C  (describing  the  various  parameters  of  the  system)  are  known  precisely. 
However,  as  mentioned  in  the  previous  section,  for  most  practical  systems,  this  is 
seldom  the  case.  Therefore,  it  will  be  assumed  that  A  G  IIIR"^",  B  G  IIIR'”^"’  and  C 
G  IIR^’*^",  where  the  notation  IIIR  stands  for  an  interval  over  the  field  of  real  numbers. 
Accordingly  A,B  and  C  are  matrices  whose  elements  (some  or  all)  are  intervals  over 
the  field  of  real  numbers.  In  a  similar  manner,  it  will  be  assumed  that  the  coefficients 
of  various  polynomials  in  the  transfer  function  (matrix)  are  intervals. 
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1.2.  Layout  of  the  Report 


The  remainder  of  this  report  is  organized  as  follows: 

1.  REVIEW  OF  INTERVAL  ARITHMETIC  TECHNIQUES:  Although  interval 
arithmetic  was  conceived  roughly  20  years  ago,  its  use  has  been  restricted  to  the 
numerical  analysis  community.  Therefore  to  make  the  project  self-contained,  a 
sufficiently  detailed  review  was  conducted. 

2.  ANALYSIS  OF  MULTIVARIABLE  SYSTEMS:  Frequently,  a  good  under¬ 
standing  and  knowledge  of  the  system  under  consideration  enables  one  to  design 
a  better  controller.  Hence  various  analysis  techniques  using  interval  arithmetic 
were  developed.  These  analysis  techniques  are  (a)  the  solution  of  systems  of  dif¬ 
ferential  equations  with  interval  coefficients  and/or  interval  initial  conditions,  (b) 
and  the  solution  of  linear  interval  equations. 

3.  DESIGN  OF  ROBUST  CONTROLLERS:  Since  the  underlying  principles  of 
interval  analysis  implicitly  account  for  a  range  of  numbers  rather  than  a  particular 
number,  it  is  easy  to  see  that  the  design,  when  it  exists,  will  automatically  be 
robust  albeit  with  some  conservatism.  In  particular,  the  problems  of  feedback 
stabilization  by  means  of  state  feedback  and  simultaneous  stabilization  for  single 
input  single  output  systems  were  addressed. 
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2.  REVIEW  OF  INTERVAL  ARITHMETIC 


2.1.  Notation 


Unless  stated  otherwise,  we  will  use  the  following  notation  in  the  rest  of  this  report 
[  a,  d  ]  interval  scalars 

d,  a  lower  and  upper  limits  of  interval  scalar 

[  A,  A  ]  interval  matrix 

A,  A  lower  and  upper  limits  of  interval  matrices 

[  a,  a  ]  interval  vector 

a,  a  lower  and  upper  limits  of  interval  vectors 

X  a  point  from  the  interval  [  x,  x  ] 

X  a  point  vector  from  the  interval  vector  [  x,  x  ] 

X  a  point  matrix  from  the  interval  matrix  [  X,  X  ] 

IR  field  of  real  numbers 

Uni  interval  scalar  assuming  value  in  IR 

niR"  interval  vector  assuming  value  in  IR" 

jjj^nxtn  interval  matrix  assuming  value  in 

[  d,  d  ]  V  [  6,  6  ]  joint  of  [  d,  d  ]  and  [6,  6  ]  =  [min{d.  6),  max{d,  6}] 


w([  A,  A  ]) 
A  0  B 
vec  (A) 


width  oi  [  A,  A]  —  A  —  A 

Kronccker  product  of  two  matrices  defined  in  eq.  (4.2) 
a  vector  of  the  columns  of  a  matrix  defined  in  cq.  (4.3) 
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2.2.  Interval  Scalars  &  Functions 


Let  n  be  the  set  of  real  compact  intervals  [a,  a  ],  a,  a  6  IR.  General  interval  operations 
can  be  defined  as  [1],  [27]: 

Definition  2.1:  Let  *  be  a  binary  operation  on  the  set  of  real  numbers  IR.  If  [  a,  a  ], 

[  6,  6  ]  6  niR,  then 

[  d,  d  ]  *  [  6,  6  ]  =  {a*b  :  a  €  [  d,  d  ],  d  €  [  fc,  6  ]},  [  d,  d  ],  [  6,  6  ]  6  I,  (2.1) 

where  *  stands  for  ,  •  and  /.  Further  the  operator  /  is  only  defined  for  the 
operations  [  d,  d  ]/[  6,  6  ],  0  ^  [  6,  6  ]. 

It  is  clear  from  the  above  definition  that 

a*b  6  [d,  d]*[fe,  S]  (2.2) 

defined  as  the  inclusion  ‘principle  of  interval  arithmetic.  It  can  be  interpreted  to  mean 
that  the  sum,  difference,  product  and  quotient  of  the  reals  (possibly  unknown)  are 
contained  in  the  sum,  difference,  product  and  quotient  of  the  including  intervals  which 
are  known  precisely. 

The  real  numbers  can  be  defined  as  point  intervals, 

a  =  [  d,  d  ]  =  [  d,  d  ]  =  [  a,  a  ].  (2.3) 

With  the  above  definition  of  point  intervals,  we  can  define  the  operations  be¬ 
tween  point  intervals  and  intervals  with  finite  width.  Further,  the  usual  arithmetic 
operations  are  carried  over  to  intervals  as  follows: 


Definition  2.2:  Let  [  d,  d  ],  [  6,  6  ]  G  HIR,  then 
(d,  d]  +  [i,  d]  =  [d  +  d,  d  +  d] 

[d,  d]  —  [d,  d]  =  [d  —  6,  d  —  d] 

.  -  .  ,  (2-4) 

[  d,  d  ]  •  [  d,  d  ]  =  [min(dd,  dd,  dd,  dd),  max(dd,  dd,  dd.  dd)] 

[d,  d]/[d,  d]  =  [d,  d]  • 

Equation  (2.4)  clearly  indicates  that  subtraction  and  division  in  interval  arith¬ 
metic  are  not  the  inverse  operations  of  addition  and  multiplication.  In  fact  this  is 
one  of  the  main  properties  of  interval  arithmetic  that  distinguishes  it  from  the  real 
arithmetic.  Some  of  the  properties  of  interval  operations  that  will  be  useful  in  the 
sequel  are: 


0  $?  [  d,  d  ]. 


Definition  2.3:  Let  [  d,  d  ],  [  d,  d  ],  [  c,  c  ]  €  HIR,  then 

[  d,  d  ]  -f  [  d,  d  ]  =  [  d,  d  ]  -1-  [  d,  d  )  (commutativity) 

[  d,  d  ]  -f-  ([  d,  d  ]  -f  [  c,  c  ])  =  ([  d,  d  ]  -f-  (  d,  d  ])  -I-  [  c,  c  ]  (associativity) 

[  d,  d  ]  •  [  5,  d  ]  =  [  d,  d  ]  •  [  d,  d  ]  (commutativity) 

[  d,  d  ]([  d,  d  ]  •  [  c,  c  ])  =  ([  d,  d  ]  •  [  d,  d  ])[  c,  c  ]  (associativity) 

[  d,  d  ]([  d,  d  ]  -f  [  c,  c  ])  C  [  d,  d  ][  d,  d  ]  -1-  [  d,  d  ][  c,  c  ]  (subdistributivity) 

a([  d,  d  ]  -I-  [  c,  c  ])  =  a[  d,  d  ]  -t-  a[  c,  c  ]  o  E  IR 

[  d,  d  ]  C  [  5,  6  ],  [  c,  c  ]  C  [  J,  d]  =>  [  d,  d  ]  *  [  c,  c  ]  C  [6,  d  ][  J,  d]  . 


(2.5) 

The  last  property  is  called  inclusion  isotonicity  of  interval  operations. 

Note  that  all  these  operations  reduce  to  real  arithmetic  when  the  intervals  are  degen¬ 
erate  intervals  or  point  intervals  i.e.,  d  =  d  and  b  =  b.  In  particular  in  the  division 
operation  with  point  intervals,  0  G  [  d,  d  ]  if  and  only  if  d  =  d  =  0. 
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Definition  2.6:  The  distance  between  two  intervals  [  d,  a  ]  and  [b,  6  ]  e  IIIR  is  defined 
as 

di[a,  d],[6,  6])  =  max{|d-6|,  |a  -  6|}.  (2.8) 

It  can  be  proved  that  the  map  d  introduces  a  metric  in  IIIR.  Note  that 

ri([  d,  d  ],  [  6,  6  ])  >  0  and  d([  d,  d  ],  [  d.  6  ])  =  0  <=>  [  d,  d  ]  =  [  6,  6  ] 
d{[  d,  d  ],  [  h,  h])  <  d([  d,  d  ],  [  r,  r  ])  +  d([  6,  b  j,  [  c.  c  ]). 

1  ] 


(2.9) 


Notice  that  (2.9)  is  the  triangle  inequality  and  may  be  verified  as  follows: 

d([  d,  d  ],  [  c,  c  ])  +  d{[  h,  6  ],  [  c,  c  ])  =  max{|n  —  c],  \a  —  c|}  +  rnax{|6  —  ci.  |/)  —  < 

>  max{|d  — cj  4-  |t  — cj,  |d  -  cj  +  \h  -c’l} 

>  max{|d  —  6j,  |d  —  b\} 

=  d([  a,  a  ],[  h  ]).  1^2.10) 

The  metric  for  real  intervals  defined  above  is  the  Hausdorff  metric  on  IIIR  and  is  a 
generalization  of  the  distance  between  two  points  in  a  metric  space.  Fnrtner,  it  has 
been  shown  that  the  metric  space  is  complete.  Existence  of  metric  in  IIR  makes  it 
a  topological  space  and  therefore  the  concepts  of  convergence  and  continuity  may  be 
used  in  the  usual  manner.  In  particular,  a  sequence  of  intervals  {[  d,  d  k  = 

0, 1, . . . ,  oo,  converges  to  an  interval  [a,  a]  if  and  only  if  the  sequence  of  the  bounds 
of  the  individual  members  of  the  sequence  converge  to  the  corresponding  bounds  d 
and  d  i.e. 

lim  [  d,  d  =  [d,  d]  <=>  (  lim  d^*'^  =  d  and  lim  d^*^  =  dV  (2.11) 

k—^oo  y<:— *oo  It— >00  j 

Theorem  2.7:  Given  a  sequence  of  intervals  {[  d,  a  assume  that  the  condition 

[  d,  d  D  [  d,  d]^’^  D  [  d,  d  D  •••  (2.12) 

is  satisfied.  Then,  the  sequence  of  intervals  {[  d,  d  converges  to  the  interval 

[d,  d]  =  n^^Jd,  d]^''^  ■ 

The  above  results  immediately  leads  to  very  useful  notions  associated  with  func¬ 
tions  in  interval  arithmetic.  They  are  stated  without  proof. 
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Theorem  2.8:  The  operations  *  G  between  intervals  axe  continuous.  B 


Theorem  2.9:  Let  /  be  a  continuous  function  and  let 


/([x,  i])  =  [  min  f(x),  max  /(x)  ], 

i€li,  r6(i,  r) 


(2.13) 


then  /([  X,  X  ])  is  a  continuous  interval  expression. 


2.3.  Interval  Matrices 

Some  of  the  results  of  point  matrices  can  be  easily  extended  to  interval  matrices. 
However,  due  to  the  nature  of  interval  operations,  several  of  them  cannot  be  extended 
directly  to  interval  matrices.  In  next  few  paragraphs,  we  discuss  some  properties  of 
interval  matrices  and  operations  with  interval  matrices  that  will  be  of  use  in  the 
sequel. 

Definition  2.10:  Two  interval  matrices  [  A,  A  ]  and  [  B,  B  ]  G  are  equal, 

i.e.,  [  A,  A  ]  =  [  B,  fi  ]  if  and  only  if  every  element  of  the  two  matrices  are  equal  t.e., 

[  dij  ]  =  [  bij,  bij  ]. 

Definition  2.11:  Given  [  A,  A  ]  and  [  B,  B  ]  G  then 

[  A,  A  ]  C  (  B,  B  ]  ^  [  aij,  d.y  ]  C  [  6.,  ],  (2.14) 

i  =  1, 2, . . . ,  n  and  j  =  1,2,  ...,m. 

The  general  operations  between  interval  matrices  can  be  defined  as: 

Definition  2.12:  ADDITION,  SUBTRACTION  -  Given  [  A,  A  ]  and  [  B,  B  ]  G  HIR"''"', 
with  elements  [  d,j,  ]  and  [bij,  bij  ]  respectively,  then 

[A,  A]  ±  [  B,  B  ]  =::  (  [  a.;,  «,;  ]  ±  [  ft,;,  ]  )  (2.15) 
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where  the  right  hand  side  denotes  the  (i,_;)-th  element  of  the  sum  or  difference  of  th( 
two  matrices  on  the  left  hand  side. 


Multiplication  -  Let  [  A,  A  ]  e  and  [  B,  B  ]  €  with  elements 

[  Oij,  Oij  ]  and  [  6,j  ,  bij  ]  respectively,  then 


[A,  A][B,  B] 


^ik  ][  brj ,  bf-j  ] 


j 


(2.16) 


defines  the  product  of  two  interval  matrices.  Interval  matrix  and  interval  vector 
product  can  be  defined  in  an  identical  manner. 

The  next  theorem  outlines  some  useful  properties  of  interval  matrices: 


Theorem  2.13;  Let  [A,  A  ],  [  B,  B  ]  and  [  C,  C  ]  be  interval  matrices  of  appropriate 
dimensions,  then 

[  A,  A  j  +  [  fi,  B  I  =  [  B,  B  I  +  [  A,  A  I  (commutativity) 

[A,  A]  +  ([B,  B]  +  [C,  C])  =  ([A,  A]  +  [B,  B])  +  [C,  C] 

(associativity) 

([  A,  A  ]  +  [  B,  B  ))[  C,  C  )  C  [  A,  A  ](  C,  C  ]  +  [  A,  A  ][  C,  C  ] 

(subdistributivity) 

[  A,  A]  +  0  =  [  A,  A]  0  is  a  null  matrix 

[A,  A]I  =  [A,  A]  I  is  an  identity  matrix 

A([B,  B]  +  [C,  C])  -  A[B,  B]  +  A[C,  C]  A  a  nominal  matrix 
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3.  SOLUTION  OF  INTERVAL  DIFFERENTIAL 
EQUATIONS 

The  solution  of  high  order  ordinary  differential  equations  (ODE)  is  a  standard  compo¬ 
nent  of  all  control  design  and  analysis  techniques  and  extremely  efficient  techniques 
exist  for  the  solution  of  the  same  [7],  [19].  However,  frequently  there  is  some  un¬ 
certainty  in  the  exact  knowledge  of  the  initial  conditions.  This  in  conjunction  with 
uncertainty  in  the  parameter  values  can  lead  to  misinterpretation  of  the  actual  solu¬ 
tions.  With  that  in  mind,  in  this  section,  we  present  the  solution  of  the  differential 
equations  with  interval  coefficients  and  interval  initial  conditions.  The  aim  of  this 
section  is  to  estimate  the  inclusion  of  the  solution  set  for  such  differential  equations. 
It  should  be  clearly  stated  here  that  the  rest  of  this  section  has  been  adopted  from  a 
paper  entitled  “An  Interval  Method  for  Systems  of  ODE”  that  appeared  in  the  Pro¬ 
ceedings  of  1985  International  Symposium  on  Interval  Mathematics.  Since  the  result 
is  an  extremely  important  one  and  not  easily  accessible,  it  has  been  rewritten  (with 
some  modifications)  and  included  in  this  report. 

Consider  the  following  system  of  n  ODE’s: 

=  f(t,  [c,  c],x(0),  x(to)  =  [xo,  Xo  ]  (3.1) 

where  c  G  [c]  =  [  c,  c  ]  G  IIIR^  and  Xq  G  [xq]  =  [  Xq,  Xq  ]  G  IHR"  are  the  system  param¬ 
eters  and  initial  conditions  respectively.  Note  that  the  uncertainty  in  the  parameters 
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as  well  as  that  in  the  initial  conditions  is  asstimed  to  lie  within  the  correspoodiag 
interval  vectors. 

Let  the  time  interval  on  which  we  wish  to  deterniine  the  solution  of  Uie 
problem  be  denoted  by  [  t,  t  ]  =  [toJf]-  Then,  the  prolilc-in  that  we  address  is  that  of 
determining  an  enclosure  [  s,  s  |  of  the  set  {x(f)}  of  all  solutions  of  (3.'' )  (assuinuig 
that  the  solutions  exist)  on  the  desired  time  interval.  Equivalently,  we  determine 
[  s,  s  ]  such  that 

s(t)  <  x(t)  <  s(t)  Vx  (3.2) 

for  each  t  €  Let  f  define  a  mapping 

f  :  f [f,  f]  X  [c,  c]  X  [a,  d],  [d,  d]eniR".  (3.3) 

where  [  a,  a  ]  is  the  domain  of  x.  We  impose  the  following  conditions  on  the  function; 

1.  f(-)  is  continuous  with  respect  to  c 

2.  Lipschitzianf  with  respect  to  x  and  t 

3.  and  isotone  with  respect  to  x. 

Then,  the  following  algorithm  will  compute  the  required  enclosures. 

Algorithm  ODE:  Compute  Enclosure  of  a  System  of  ODE 

Let  /i  >  0  be  a  sufficiently  small  step  such  that  —  tg  +  kh  ^  [i^  f  ].  where  k  is 
an  integer  and  let  A:  =  0, 1, ,  k. 

Let  s(fo)  =  Xo  and  s(fo)  =  Xo; 

Let  F,([f],[.D])  =  {/,(f,a-)  :  t  €  (T],x  €  p]},  F.df],  [D])  -  {/.(t.,r)  :  t  e 
[T],x  G  [.D]},  i  =  1, 2, . . . ,  n  for  each  [T]  C  [  t,  f  ]  and  [!)]  C  [  a,  a  ]. 
t  A  function  satisfying  Lipscliitz  condition:  \f{x,y)  —  f{x.y)\  <  (\.t  —  .f|,  is  said  to 
be  Lischitzizm. 
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Assume  further  that  we  have  already  computed  s(tk)  and  s(tfc),  such  that  Si{tk) 
<  Xi(tfc)  <  i  =  1,2,  Then  the  following  recursion  computes  the 

vectors 

s(t)  =  [ii  S2  in]^  and 
s(<)  =  [il  52  •  •  •  -Sn]^, 
in  the  interval 
upper  bound  s{t): 

=  [d.,d,],  i  =  1,2,  ...,n 
for  r  =  0, 1, . . .  ,r, 

i  =  l,2,...,n 

=  i,(4)  V  (si(tk)  +  V  (si{tk)  +  i  =  1,2,  ...,n 

Si(t)  =  Siitk)  +  (t -tk)qi^\  i  =  l,2,...,n,  t  €  [<fc,tjb+i] 


lower  bound  s(t): 

=  [di,d,],  i  =  1,2,. . .  ,n 
for  r  =  0, 1, . . . ,  f, 

j  =  1,2,. ..,n 

=Si{tk)y  (^Si(tk)  +  V  (.8i(4)  +  *  =  l,2,...,n 

Si{t)  =  Siitk)  +  (t -tk)p\’'\  i  =  l,2,...,n,  <  €  [4,4+1] 


Effectively,  the  algorithm  is  computing  the  following  quantities: 

1.  the  intervals  /,(t,  [  c,  c],x(t))  =  {/,(f,  c,  x(f )),  c  €  [  c,  c]},  i  =  1,2,. ..,n., 

for  each  value  of  f  €  ^nd  the  vector  x(t)  €  [  d,  d  ].  The  functions 
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and  F;  ^4,  4’  . '  =  1'2 .  ..  v 

spectively,  denote  the  end  points  of  each  of  the  interval. 

2.  and  the  intervals 

Fi([f],[D])  =  :  t  e  [r],x  6  [D]  and 

F,([T],[D])=  {/.(t,x):te[f],xep], 

i  =  1, 2, . . . ,  n  for  each  [T]  €  [  f,  t  ]  and  [f)]  e  [  d,  d  ]. 

At  each  time  step,  the  solution  interval  is  obtained  by  determining  the  smallest 
interval  that  will  include  the  largest  deviations  in  the  values  of  x,  G  [  d,.  d,  |.  It  is 
should  be  pointed  out  that  the  resulting  solution  vector  will  overbouiid  the  soluMons 
generated  by  all  possible  variations  in  [  t,  t  ],  [  c,  c  ]  and  [  d,  d  ]. 

To  prove  the  convergence  of  the  above  algorithm,  we  proceed  as  follows: 

Theorem  3.1:  For  any  non-negative  integer  r. 


,(r+l)  .(r)  j  -(r+1)  -(r) 

C  2,  and  z-  C  2, 


(3.4) 


Proof:  See  Appendix  3. A. 

Next,  it  will  be  shown  that  s  and  s  are  indeed  bounds  for  the  solution  set.  To 
show  this  we  will  establish  that  the  sohitions  x(t)  are  bounded  as  follows: 

Theorem  3.2:  Given  the  solution  set  {x(<)},  it  is  bounded  as 


s(0  <  {x(t)}  <  s(<),  tk  G  [tkJk+^].  (3.5) 


Proof:  See  Appendix  3.B. 

Knowing  the  lower  and  upper  bounds  of  the  time  response  of  a  linear  interv^al 
differentied  equation,  feedback  can  be  applied  to  improve  the  transient  and/or  steady 
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state  response  of  the  system.  Study  of  feedback  vis  a  vis  improvement  of  time  response 
is  of  considerable  importance  and  is  a  possible  direction  in  which  to  extend  the  results 
presented  in  this  report.  It  is  conjectured  that  similar  results  can  be  derived  for 
bounding  the  gain  and  phase  plots  for  rational  functions  with  interval  coefficients. 

Appendix  s.a. 

Proof  of  Theorem  3.1:  The  proof  is  by  induction.  For  r  =  0,  we  have 

=Si{tk)y  (siitk)+p\°^hjy  (^Si(tk)  +  q\°^h^,  i  =  1,2,. ..,n.  (3.6) 

Since,  3i(tk)  6  i  =  1,2, ...  ,n,  we  can  make  the  increment  h  sufficiently  small 


such  that  C  [di,di]  = 

Next,  assume  that  C  for  some  r  >  2. 

Then,  in  the  interval  [<*,<*+1], 

since  Fi  is  an  inclusion  isotone, 

pi"'. 9!’'’] 

,  i  =  l,2,...,n. 

(3.7) 

Therefore, 


p'r"  <  pV  <  si”’  <  9!'"” 


and, 


Equivalently, 


Si{^k)^pY  <  s,{tk)  +  p\''h 

Siitk)  +  >  Si{tk)  + 

^  •‘i  • 


(3.8) 


(3.9) 


(3.10) 
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In  a  similar  manner,  it  can  be  shown  that 


.(r+l) 

• 


Appendix  s.b. 


Proof  of  Theorem  3.2:  For  any  non- negative  integer  r, 


r(0  =  ^  t  €  (4, <*4.,],  Sj  6  Zj'"',  j  =  1.2,..., 


n. 


13.11) 


From  Theorem  3.1,  it  is  known  that  C  then  for  each  t  £ 


.(.--H)  ^  -(r) 


(3.12) 


Sj{t)  =  Sj(tk)  +  pY\*-*k)  €  C  .  J  =  1,2 . n.  (3.13) 

Therefore, 

^i(0  <  /.(t,3i(t),...,.s„(0),  *  =  l,2,...,n.  (3.14) 

Following  the  proof  in  an  analogous  manner,  it  can  be  established  that 

(<),.■■  ,'Sn(t)),  i  -  1,2 - ,n.  (3.15) 


If  we  assume  that  x(t)  is  an  arbitrary  solution  of  the  interval  ODE,  corresponding 
to  some  c  £  [c]  emd  some  Xo  €  [Xo],  then, 

s(t)  <  f(t,s(t))  <  f(t,[c],s(0) 

lit)  >  f(^,s(<))  >  f(t,[c],s(t))  (3.16) 

s(<o)  =  Xo  <  Xo  <  Xo  =  S(t). 

The  above  relations  can  be  summarized  as 

s{t)  <  f(t,[c],s(0), 
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x(i)  =  f(<,[c],x(0), 

s{t)  >  f(t,[c],s(i)), 

sf^o)  <  Xo, 

x{to)  =  Xo  and 
s(<o)  >  Xo. 

Under  the  assumption  that  f  is  a  quasi-isotone  in  x, 

s(t)  <  {x(<)}  <  s{t). 
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4.  A  SUFFICIENT  CONDITION  FOR  STABILITY  OF  A 
CLASS  OF  INTERVAL  MATRICES 

In  this  Section  we  derive  a  sufficient  condition  for  ensuring  the  stability  of  a  class  of 
interval  matrices  called  -matrices  and  study  its  applications.  The  Lyapunov  crite¬ 
rion  of  stability  of  a  matrix  is  extended  to  interval  Ad-matrix  matrix  case.  However, 
in  extending  to  the  interval  case,  seme  conservatism  invariably  enters  the  analysis 
thereby  rendering  the  criterion  only  sufficient.  Some  applications  of  the  stability  of 
Ad -matrices  to  dynamical  systems  are  presented. 

It  is  well  known  that  for  systems  described  by  nominal  (intervals  with  zero  width) 
parameters, 

Theorem  4.1;  Let  V(x)  =  x^Px,  where  P  is  a  solution  to  the  matrix  Lyapunov 
equation 

A^P  -b  PA  =  -Q.  (4.1) 

Then  a  real  matrix  A  is  a  stability  matrix  if  and  only  if,  for  any  given  real  symmetric 
positive  definite  matrix  Q,  the  solution  P  is  also  symmetric  positive  definite.  ■ 

The  solution  of  a  matrix  Lyapunov  equation  can  be  obtained  using  existing  nu¬ 
merically  reliable  software  for  nominal  systems.  However,  these  techniques  cannot 
be  extended  to  intervsil  matrices.  Fortunately,  we  can  use  the  following  alternate 
technique  for  solving  the  matrix  Lyapunov  equation. 
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Define  the  Kronecker  product  of  two  matrices  A  =  (aij)  €  and  B  €  IR**  ® 


as  A  (2)  B  6  IR"*'’^"®,  where  the  product  is  given  by 


A®B 


UiiB  012^ 

021 B  022  S 


oi„B 

02nB 


LdmlB  flmnBJ 


Further,  we  can  define 


vec  (P) 


Pi  ' 
P2 


Lp„J 


where  pi  is  the  i-th  column  of  P.  On  noting  that 


vec  (ABC)  =  (g)  A  vec  (B) 


(4.2) 


(4.3) 


(4.4) 


and 

A^P  +  PA  =  A^PI  +  IPA  =  -  Q,  (4.5) 

t2Lking  vec  (•)  on  both  sides,  we  get 

(l  0  A^  +  A^  0  l)  vec  (P)  =  —  vec  (Q).  (4.6) 

Equation  (4.6)  can  be  solved  by  Gaussian  elimination  provided  (I  0A^  +  A^0I) 
is  nonsingular.  By  the  properties  of  Kronecker  sums,!  it  is  clear  that  the  matrix  on 
the  left  hand  side  will  be  nonsingular  if  and  only  if  A,  +  Aj  ^  0  for  all  eigenvalues  of 
A.  Clearly,  if  A  is  stable  (all  eigenvalues  in  left  half  plane),  then  the  above  condition 
is  satisfied  and  a  unique  solution  exists. 

t  The  eigenvalues  of  A0B  are  the  mn  numbers  A^/ij,  i  =  1, 2, . . . ,  n  j  =  1, 2, . . . ,  m, 
where  Ai  and  fij  are  respectively  the  eigenvalues  of  A  and  B.  Further,  the  eigenvalues 
of  A  0  I„  -f  Im  0)  B  (called  the  Kronecker  sum)  are  the  numbers  Xi  +  fij. 
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For  the  case  of  interval  systems,  we  cannot  solve  for  P  directly  because  it  is  not 
easy  to  compute  the  inverse  of  an  interval  matrix.  The  best  that  can  be  done  is  obtain 
an  enclosure  of  the  solution  vector  -  vec  (P).  We  next  investigate  how  to  compute 
these  enclosures. 

4.1.  Solution  of  Linear  Interval  Equations 

In  this  Section  we  study  the  problem  of  solving  linear  interval  equations.  These 
will  play  a  significant  role  in  developing  the  Lyapunov  stability  criterion  for  interval 
Af-matrices.  In  the  literature,  a  linear  interval  equation  with  a  coefficient  matrix 
[  A,  A  ]  and  a  right  hand  side  [  b,  b  ]  is  defined  as  a  family  of  the  linear  equations 
characterized  by 

Ax  =  b  (4.7) 

where  A  €  [  A,  A  ]  and  b  €  [  b,  b  ].  Note  that  (4.7)  implicitly  assumes  that  A  exists 
for  all  A  €  [  A,  A  ].  The  solution  set  of  (4.7)  may  be  defined  as 

S  A,  A  ],  [  b,  b  :=  |x  6  IR”| Ax  =  b  for  some  A  €  [  A,  A  ],  b  G  [  b,  b  ]  |  . 

(4.8) 

An  obvious  way  to  compute  such  a  vector  [  x,  x  ]  is  to  extend  Gaussian  elim¬ 
ination  to  systems  of  linear  equations  with  interval  coefficients.  For  simplicity  of 
presentation,  assume  that  the  (interval)  elements  of  [  A,  A  ]  are  defined  as  a,j  and 
the  elements  of  [  b,  b  ]  are  defined  as  6;.  Then  starting  with  the  following  coefficient 
tableau: 


an 

0-12 

n 

a>i  022 

^2n 

^nl 

«n2 

n 

bn 
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2ind  employing  the  relations: 
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for  which 


Theorem  4.3:  Let  1  <  n  <  2  and  assume  that  the  n  x  n  interval  matrix  [A,  A  ]  = 
([  dij,  djj  ])  does  not  contain  a  singular  matrix  A,  then  the  Gaussian  algorithm  can 
always  be  carried  out. 

Proof:  See  Appendix  4.B. 

Unfortunately,  the  above  result  is  true  for  n  x  n  interval  matrices  where  1  <  n  <  2 
only  and  cannot  be  extended  to  the  case  when  n  >  3.  Note  that  Gaussian  elimination 
can  be  performed  for  the  several  special  classes  of  matrices.  One  of  the  more  important 
ones  being  At-matrices.  In  the  next  section  we  will  show  that  if  the  interval  matrix 
is  also  an  Al-matrix,  then  the  solution  of  linear  interval  equations  is  considerably 
simplified. 


26 


4.2.  Convex  Hull  of  m-Matrices 


We  will  consider  the  case  of  square  interval  matrices.  The  results  can  also  be  extended 
to  non-square  interval  matrices,  with  applications  in  solving  interval  least  square 
problems. 

Definition  4.4:  A  square  matrix  A  is  called  an  Af-matrix  if  a,j  <  0  for  i  ^  j  and  Au 
>  0  for  some  positive  vector  u  €  IR”. 

An  important  property  of  nonsingulax  Af-matrices  is  that  every  element  of  their 
inverse  is  non-negative. 

As  mentioned  earlier,  the  computation  of  [  A,  A  ]  for  a  regularf  interval  ma- 
trix  is  a  difficult  problem.  The  known  methods  have  a  computational  complexity 
exponential  in  n  (the  order  of  the  matrix).  But  for  the  important  special  case  of 
Af-matrices,  we  have  the  following  explicit  result  [28],  [29]: 

Theorem  4.5:  For  an  n  X  n  interval  matrix  [A,  A  ],  if  the  matrices  A  and  A  are 
regular  and  A~^  >  0,  A~^  >  0,  then  [  A,  A  ]  is  regular  and 


[A,  A]  ‘  =  [A-^  A-' 


(4.16) 


Proof:  See  Appendix  4.C. 

4.3.  Stability  of  Interval  at-Matrices 

Researchers  have  attempted  to  derive  necessary  and  sufficient  conditions  for  estab¬ 
lishing  the  Hurwitz  or  Schur  stability  of  interval  matrices.  However,  invariably,  the 
results  have  been  shown  true  for  only  a  smaller  class  of  matrices  or  the  conditions 
t  If  all  matrices  A  G  [  A,  A  ]  have  full  rank,  then  [  A,  A  ]  is  a  regular  matrix 
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have  been  proved  only  to  be  sufficient.  Notably,  Bialas  [10]  presented  a  result  that 
derived  a  necessary  and  sufficient  condition  of  stability  of  general  interval  matrices. 
This  result  was  later  disproved  by  several  researchers  e.g.,  Barmish  and  Hollot  [6]  and 
Kaxl  et  al.,  [20].  A  sufficient  condition  for  stability  of  interval  matrices  is  derived  by 
Heinen  [17]  and  Yedavalli  [43].  In  [15],  the  authors  obtain  the  conditions  for  stability 
of  a  special  class  of  matrices.  Other  related  work  has  been  in  obtaining  the  perturba¬ 
tion  bounds  on  the  elements  of  a  nominal  stable  matrix.  These  bounds  can  be  found 
in  the  work  of  Patel,  Toda  and  Sridhar  [31],  [32]  and  Yedavalli  [42]. 

In  this  Section,  we  will  derive  the  condition  for  the  stability  of  interval  ,M- 
matrices.  The  following  important  properties  of  the  A4-matrices  will  be  used  in  the 
sequel; 

Theorem  4.6:  POSITIVE  STABILITY  [8]  A  is  a  nonsingular  Al-matrix  if  and  only  if 
there  exists  a  solution  W  to  the  matrix  equation 

AW  +  WA^  -  Q  (4.17) 

where  W  and  Q  are  symmetric  positive  definite  matrices.  B 

Theorem  4.7:  POSITIVITY  OF  PRINCIPAL  MINORS  [8]  A  is  a  nonsingular  .M-rnatrix 
if  and  only  if  A  is  nonsingular  and  x  0;  y  =  Ax,  then  for  some  subscript  .r,  7^ 
0,  and  Xiyi  >0.  B 

Theorem  4.8:  If  A  is  a  nonsingular  Al-matrix  then  (I  ®  A^  -|-  A^  (8)  I)  is  also  a 
nonsingular  Af-matrix. 

Proof:  See  Appendix  4.D. 
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The  above  result  is  fairly  important  for  the  sequel,  hence  it  is  illustrated  by  means 


of  a  simple  numerical  example. 


Example  4.1:  Consider  the  (2  x  2)  non-singular  -matrix 


A  = 


Forming  the  Kroneckar  products,  we  get 


■  5 

-2 

0 

0  - 

■  5 

0 

-2 

0  ■ 

-1 

4 

0 

0 

,  A^  ®  I  = 

0 

5 

0 

-2 

0 

0 

5 

-2 

-1 

0 

4 

0 

.  0 

0 

-1 

4  . 

.  0 

-1 

0 

4  . 

Therefore, 


I(8)A^  -I-  A^(8)I 


The  inverse  of  the  above  matrix  is 


■10-2-2  0  ■ 
-19  0-2 

-10  9-2 

.0-1-1  8  . 


-1.0494D-01 
1.23461?  -  02 
1.2346I>  -  02 
.3.0864D-03 


2.4691I>  -  02 
1.1728I>  -  01 
6.1728I1-03 
1.54321?  -  02 


2.4691Z?  -  02 
6.1728Z?  -  03 
1.1728i?  -  01 
1.54322?  -  02 


1.23462? -02- 
3.08642?  -  02 
3.08642?  -  02 
1.32722? -01. 


Thereby,  verifying  that  I  (gi  A^  -I-  A^  0  I  is  a  non-singulax  A4-matrix. 


Next,  we  determine  the  condition  for  the  positive  stability  of  interval  Al-matrices. 
By  definition,  all  eigenvalues  of  [A,  A]  are  positive  (i.e.  the  matrix  is  positive 
stable).  Therefore,  to  obtain  a  Lyapunov  function  to  establish  positive  stability  of 
interval  matrices,  we  need  to  find  an  interval  positive  definite  matrix  as  a  solution  of 
the  following  Lyapunov  equation 


[  A,  A  ]'"’[  P,  P  ]  4-  [  P,  P  ][  A,  A  ]  =  Q.  (4.18) 
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where  Q  is  positive  definite  matrix.  Then,  if  each  P  €  [  P.  P  ]  is  a  jiositive  definite 
matrix,  [A,  A  ]  is  a  positive  stable  matrix. 

Equivalently,  assuming  that  [A,  A  ]  is  an  interval  ,,V(-matrix,  we  form  the  Kro- 
necker  sum  given  by  (4.2)  for  internal  matrix  case  and  solve  the  following  linear  interval 
equation  for  vec  (P): 

(I®[A,  A]^  +  [A,  A]^®I)  vec{[P,  P])  =  vec(Q).  (4.19) 
Denoting  the  v?  x  matrix  on  the  left  hand  side  by  (  A.  A  ].  we  have 

[  A,  A]  vec  ([  P,  P  ])  =  vec  (Q).  (4.20) 

where  vec  (Q)  >  0.  Now,  using  the  properties  of  interval  Al-matrices, 

vec  ([  P,  P  ])  =  [  v4,  ,4  ]  ^vec  (Q) 

=  [  A~^  ]vec  (Q) 

=  [  A“*vec  (Q),  A“'vec(Q)].  (4.21) 


The  interval  matrix  [  P,  P  ]  can  now  be  easily  constructed.  Clearly,  [  P,  P  ]  is  a 
symmetric  matrix. 

To  verify  whether  it  is  positive  definite,  we  use  the  result  developed  by  Shi  and 
Gao  [38].  Following  the  notation  used  in  [38],  let  U  =  (u,j)  and  V  =  (t’,j  )  G 
be  symmetric  matrices  and  <  Vij,  i,j  =  1,2, .. .  ,n.  Denoting 

£"(U,V]  =  {A  =  (a,^)  €  <  n.j  <  a,j  =  =  1,2, - n} 

A:"[U,V]  =  {S  =  (s,j)  €  =  7i,^-  or  Sij  =  u.j;  ?,  j  =  1,2 . n). 

(4.22) 
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Note  that  while  the  first  equation  represents  a  set  of  matrices  with  infinite  members, 
the  second  equation  denotes  a  set  of  edge  matrices.  According  to  [38],  an  interval 
matrix  £"[U,  V]  (/C”[U,V])  is  termed  a  positive  definite  matrix  if  A  (S)  is  positive 
definite  for  every  A  €  £”[U,V]  (S  €  A:"[U,V]).  Based  on  the  above  notation,  the 
following  result  wcis  established  in  [38]: 

Theorem  4.9:  The  set  £"[U,  V]  consists  of  positive  definite  matrices  if  and  only  if 
the  set  ^"[U,  V])  consists  of  only  positive  definite  matrices.  B 

Using  the  above  result,  now  it  is  clear  that  [  P,  P  ],  the  solution  of  interval  Lya¬ 
punov  equation  is  positive  definite  if  and  only  if  all  the  vertex  matrices  are  positive 
definite.  To  establish  the  positive  definiteness  one  would  have  to  verify  positive  defi¬ 
niteness  of  all  vertex  matrices.  Note,  however,  the  matrix  A  is  (4.20)  is  over-bounded, 
hence  positive  definiteness  of  [  P,  P  ]  is  only  sufficient  for  positive  stability  of  [  A,  A  ]. 
Note  further  that  if  the  entire  interval  Lyapunov  equation  is  premultiplied  with  nega¬ 
tive  identity,  parallel  results  can  be  obtained  for  Hurwitz  stability  of  inverse  negative 
matrices  (matrices  whose  inverse  has  all  negative  elements). 

4.4.  Applications  of  Stability  of  ai-Matrices 

In  this  Section,  we  study  two  applications  of  the  stability  conditions  developed  above. 
The  first  application  is  in  the  study  of  robust  stability  of  composite  dynamical  sys¬ 
tems  and  the  second  application  is  in  the  study  of  robustness  of  dynamical  Leontief 
economic  models. 

4.4.1.  ROBUST  STABILITY  OF  DECENTRALIZED  SYSTEMS 

The  earliest  application  of  ^Vf-matrix  theory  to  composite  dynamical  system  stability 
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problem  was  studied  by  Gurjic  and  Siljak  [IS);  and  Araki  [2].  In  the  sequel  we  will 
follow  the  notation  used  by  Gurjic  and  Siljak.  Consider  a  continuous  time  system 
described  by  a  set  of  n  vector  differential  equations  [37]: 


/IX 

=  g.(x.,t)  +  hCA’,^,  i  =  1,2,. ...k  (4.23) 

where  x,  €  K”'  ,  gi  :  IR”‘  x  IR  — >  IR"’  ,  hi  :  IR"*  x  IR"*  x  •  ■  •  x  IR"*  x  IR  IR"'  and 
.X  =  {xi,...,Xit}. 

Assume  that 


g,(0,t)  =  0  V  t  €  IR  (4.24a) 

h.(0,...,0,t)  =  0  V  t  e  IR  (4.246) 

i.e.,  the  null  state  is  the  equilibrium  state.  Clearly,  the  above  system  can  l:)e  regarded 
as  a  composite  system  consisting  of  k  subsystems  give  by 

c/x ' 

=  g,(x,,t)  +  y,.  i  =  l,2,...,k  (4.25) 

where  yi  denotes  the  interconnecting  relations  between  the  various  stations.  When 
the  interaction  hi(A’,t)  =  0,  the  unforced  subsystem  is  defined  as 

^  =  g.(x,,t).  (4.26) 

In  the  subsequent  paragraphs  we  will  assume  that  the  stability  of  each  unforced 
system  has  been  established  by  determining  scalar  Lyapunov  functions  v,  :  IR  x  IR"' 
— >  IR  satisfying 


(f>ti  (||x,j|)  <  c.(/,x,)  <<;j,2{||x,| 
de,(t,x,) 


(It 


<«!>..■)  (||x, I 


and  (4.27a) 

V  t  6  IR  and  V  x,  €  IR"’.  (4.276) 


The  total  time  derivative  of  is  given  by 

Vi  =  ^Vi  +  (grad  (4.28) 

Further,  it  is  assumed  that  there  exist  bounded  functions  ^ij  :  IR  x  IR”  — ♦  IR, 

sup  l^ijl  =  Vij  <  oo  V  =  1,2,  (4.29) 

IRxlR” 

subject  to  the  following  inequalities: 

k 

(graduO'^h,  <  i  =  1,2,. (4.30) 

j=i 

Next,  we  study  the  composite  stability  of  the  system,  knowing  the  stability  prop¬ 
erties  of  individual  subsystems.  The  total  time  derivative  of  the  function  v,(A’,t) 
along  the  solutions  of  unforced  individual  subsystems  is  given  by 

iii  =  Vi  +  (gradvj)^h,.  (4.31) 

Using  the  bounds  in  (4.27)  and  (4.30),  we  get 

k 

ii  <  <^.3(|lx.l|)  +  »  =  1,2,. ..,fc.  (4.32) 

Denoting  the  vector  Lyapunov  function  Jis  v  =  [vi,  and  the  comparison 

vector  function  ^  =  [<^13,  •  •  •  ,<^<:3]^  related  as 

V  <  A$,  (4.33) 

the  elements  of  A  are  given  by  (4.32)  as 


Qtj  —  Sij  -t-  Ckij 


(4.34) 
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where,  6ij  is  the  Kronecker  symbol  and 


ai>  =  ^ij  sup  -h  max{0,  sup  (4.35) 

IR"  xIR  IR"  xIR 

Clearly,  the  stability  of  the  composite  system  is  gc  vernecl  by  the  stability  cjf  the 
dynamical  system  defined  by  (4.33).  Notice  that  in  (4.33),  A  is  an  ^Vl-nic-irix. 

Now,  to  illustrate  the  use  of  the  theory  developed  in  the  previous  Section,  let  us 
consider  the  example  in  [18]. 

Example  4.2:  Assume  that  the  system  is  composed  of  two  subsystems  defined  by 


=  g,(x,)  +  (4.3Ga) 

~  =  gi(xi,<)  4-  hi(A',t)  (4.366) 


where  A  =  [xf  x|’]^,  Xi  =  [x,,  2^12]^,  X2  =  [x2i  122]^, 


gi 

hi 

g2 

h2 


3^12 

-2xu  -  X^iXi2 

■  |[  d,  a  ]xiixf2  +  [/?,/?  ]{xl^  +  X22) 

.  [  d,  d  ]x?ixi2  +  [/?,/?  ](x|i  +  x^2)  . 
X22  cos(t)  —  X2i(x2i  +  X22)(-i  +  sin(t)) 
-X21  cos(t)  -  X22(x2i  +  X22)(2  +  sin(t)) 

■[6,  6]x2i(  X2J  +  X22  )  +  ^(Xi,t)l 

.[  6,  6  ]x22(x|,  +  X^2)  +  ^{Xl,0. 


(4.37a) 

(4.376) 

(4.37c) 

(4.37(/) 


where  |d(xi,t)|  <  [7,  7]<^«3(||xi||)  and  (  d,  d  ],  [  ,3,  /I  ].  [  7,  7  ]  and  [  6,  6  ]  are  all 
positive  interval  parameters.  The  interval  parameters  are  a  consequence  of  uncertain¬ 
ties  in  the  interconnection  mo<lels.  Selecting 


Ui(xi)  =  [0.5(2xJ,  -|-xJ2)]‘''^  V2(X2)  =  IIX2II  (4.38) 
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and  following  equations  (4.30)-(4.35),  we  get 


^1  <  “  (1  —  8[  d  ])(^13  (llxill)  +  8[ /3,  ^  ](;/)23  (IIX2II)  (4.39a) 

V2  <  2(7,  7  ])<^i3(l|xi||)  -  ^  ])<^23(l|X2ll).  (4.396) 


Equivalently, 


-1  +  8[  d,  a]  8[  /3 

2[  7,  7  ]  -1  +  [  6,  6  ] 


(4.40) 


Notice  that,  since  [A,  A  ]  in  (4.40)  is  only  a  2  x  2  matrix,  the  stability  is  guaranteed 
provided  both  diagonal  elements  axe  negative.  However,  in  general  it  will  be  difficult 
to  ensure  the  stability  conditions. 

Fortunately,  [A,  A  ]  is  an  interval  A<-matrix,  thereby,  simplifying  the  analysis. 
To  apply  interval  matrix  analysis  developed  in  previous  Section,  we  study  the  positive 
stability  of  — [  A,  A  ].  If  —[A,  A  ]  is  positive  stable,  then  [  A,  A  ]  will  have  Hurwitz 
stability. 


4.4.2.  DYNAMIC  LEONTIEF  MODELS 

One  of  the  simpler  static  economical  models  [41]  has  the  following  characteristics: 

1.  no  joint  production;  each  sector  or  industry  produces  one  and  only  one  com¬ 
modity,  i.e.,  there  is  one  to  one  correspondence  between  the  s'-ctors  and  the 
commodities 

2.  each  sector  has  only  one  technique  of  production;  to  produce  one  unit  of  com¬ 
modity  j,  sector  j  requires  Oij  units  of  commodity  i;  i,j  =  1, 2, . . . ,  n 

3.  there  are  no  production  lags 

4.  there  is  no  government  activity  or  foreign  trade 
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5.  the  model  allows  for  only  those  commodities  which  cease  to  exist  in  the  industry, 
once  they  are  used  in  the  production  process. 

With  these  assumptions  on  the  activity,  and  denoting  the  following  for  i.  j  -= 
1,2, ...  ,n: 

Xi  =  gross  output  of  sector  i 

Xij  =  amount  of  commodity  i  used  by  the  sector  j 

=  XijIXj 

Ci  —  final  demand  for  commodity  i 
the  balance  equations  for  each  sector  can  be  written  as 


n 


X.  =  +  c<,  =  1,2,.. 

n 

.  ,n 

(4.41a) 

Xi  =  +  c„ 

j=i 

(4.416) 

^  X  =  Ax  +  c, 

(4.41c) 

where  Oij  is  the  (i,  j)-th  element  of  the  matrix  A.  The  model  in  (4.41)  is  the  static 
economic  model.  The  properties  of  the  static  model  can  be  studied  from  the  above 
system  of  linear  equations. 

This  model  can  be  further  generalized  to  the  dynamic  model.  For  the  dynamic 
model,  we  modify  assumption  (2)  as  follows  —  to  produce  one  unit  of  commodity  j , 
sector  j  requires  Oij  units  as  the  current  input  of  commodity  i  and  6,y  units  as  the 
capital  input.  Then  defining 

Sij  —  the  stock  of  commodity  i  held  by  sector  j 

bij  =  Sij/ij, 
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the  continuous  time  balance  equations  for  each  sector  can  be  rewritten  as 


=  +  Cj,  i  =  1,2,  ...,n 


X  = 


n  n 

^  ^  O-ijXj  +  ^  ^  ^ij  +  Cj,  I  =  1,2,  ...,n 

>=1  i=i 

Ax  +  B—  +  c, 
at 


(4.42a) 

(4.42a) 

(4.42c) 


The  above  model  can  be  expressed  into  the  standard  state  space  for  as 

/7y 

^  =  B-\I  -  A)x  +  B-'c  (4.43) 

at 

In  general,  the  capital  input  is  known  very  accurately.  Allowing  perturbations  in 
the  input  of  the  current  commodity  ajj,  we  have  the  following  interval  plant: 

^  =  B-^(I  -  [A,  A])x  +  B-*c  (4.44) 

at 

where  I  —  [  A,  A  ]  is  a  non-singular  interval  Af-matrix,  therefore  we  can  apply  the 
results  of  Sections  4. 1-4.3  and  study  the  stability  of  dynamic  Leontief  models. 

Next,  we  present  an  example  to  illustrate  the  robust  stability  of  uncertain  Leontief 
models. 


Example  4.3:  Following  the  above  notation,  assume  that  the  current  commodity 
has  some  uncertainty  and  the  capital  input  is  known  precisely.  The  corresponding 

matrices  [  A,  A  ]  and  B  are  given  by: 

■[0.05,0.10]  [0.10,0.30]  [0.20,0.25]]  Fl  0  O' 

[A,  A]:=  [0.00,0.10]  [0.15,0.20]  [0.20,0.30]  ,  B=  0  1  0 

[[0.10,0.20]  [0.10,0.15]  [0.05,0.10]J  [o  0  1 

Then  the  system  matrix  is  given  by 

r  [-hO.OO, -t-0.95]  [-0.30,-0.10]  [-0.25,-0.20]' 

B-MI-[A,  A])  =  [-0.10,-^0.001  (-^0.80, -hO.85]  [-0.30,-0.20]  . 

[[-0.20,-0.10]  [-0.15,-0.10]  [+0.90,-F0.95] 
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"10  0“ 

Solving  the  associated  interval  Lyapunov  equations  with  Q  =  0  10,  \vv 

[O  0  1 

have 

■[0.5521,0.6988]  [0.0553,0.2432]  [0.0946,0.2238]' 

[  P,  P  ]  =  [0.0553,  0.2432]  [0.6112,  0.7380]  [0.0978,  0.2202] 

[[0.0946,0.2238]  [0.0978,0.2202]  [0.5466,0.6420] 

which  is  a  symmetric  matrix  (rounded  to  4  decimal  places).  Notice  that  it  is  diagonally 

dominant;  therefore,  one  can  expect  favorable  results.  However,  using  the  result  in 

[38],  it  was  verified  that  the  vertex  matrices  of  [  P,  P  ]  are  indeed  positive  definite. 

Therefore  the  system  exhibits  positive  stability.  Since  [  P,  P  ]  is  a  symmetric  matrix 

of  order  3,  there  will  be  2®  edge  matrices.  In  order  to  save  space,  we  have  given  below 

the  plots  of  the  eigenvalues  of  [  P,  P  ].  FIGURE  4.1  plots  the  three  eigenvalues  as  the 

elements  of  the  matrix  [  P,  P  ]  are  varied.  The  matrices  are  in  the  following  sequence: 

the  first  matrix  Vi  is  the  matrix  of  lower  bounds  of  each  interval  i.e.  flP]  =  P.  In 

the  second  matrix,  the  (1,1)  element  is  the  upper  bound  of  [  pn,  pn  ]  and  the  rest 

of  the  elements  axe  the  lower  bounds, 

'0.6988  0.0553  0.0946' 

V2  =  0.0553  0.6112  0.0978  . 

0.0946  0.0978  0.5466 

The  third  matrix  has  its  (1,2)  and  (2, 1)  elements  as  the  upper  bound  of  [  pi2<  P\2  ] 

and  the  rest  of  the  elements  are  the  lower  bounds. 

'0.5521  0.2432  0.0946' 

V3  =  0.2432  0.6112  0.0978  . 

0.0946  0.0978  0.5466 

In  this  manner  elements  are  varied  till  the  matrix  becomes  P.  All  eigenvahu's  of 
each  matrix  are  positive,  thereby  establishing  positive  definiteness  of  the  symmetric 
interval  matrix  [  P,  P  ]. 
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Figure  4.1;  EIGENVALUES  OF  VARIOUS  EDGE  MATRICES 
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Appendix  4.a. 


Proof  of  Theorem  4.2:  Assume  that  A  6  [  A.  A  ]  and  b  G  [  b,  b  ]  and  consider 
the  system  of  linear  equations 

Ax  =  b.  (4.45) 

Let  A  =  (ajj),  the  vector  b  =  (6,),  A'  =  (a'j)  and  the  vector  b'  =  (6'^),  where. 

a\j  =  aij  1  <  j  <  n 
a'ij  =  a.j  -  aij(a,i/aii)  2  <  i,  j  <  n 

b\  =  6, 

b'-  =  bi  —  b\{aii  fail)  ’2,  <  i  <  n 

a'jj  =  0  2  <  i  <  n 

where  0  ^  on.  Clearly  A'y  =  b'  has  the  same  solution  as  Ax  =  b.  Then  by  the 
inclusion  monotonicity  (DEFINITION  2.5),  of  interval  operations,  it  follows  that 

A'  €  [  A',  A'  ]  and  b'  G  [  b',  b'  ].  (4.46) 

which  proves  the  assertion.  ■ 

Appendix  4.b. 

Proof  of  Theorem  4.3:  For  n  =  1,  the  assumption  that  0  ^  [dj],  dn  ]  directly 
proves  the  assertion.  Fcjr  v  —  2,  at  h'ast  one  of  tin'  intervals  [  dj  i ,  rin  ]  and  [  a2\ ,  ] 

must  not  contain  zero.  If  this  was  not  the  case  then  there  would  exist  a  singular  matrix 
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A  €  [A,  A  ],  contradicting  the  assumption  of  the  theorem.  Assume  without  loss  of 
generality,  0  ^  [  on,  an  ]•  Using  the  Gaussian  elimination,  we  have 


022  =  a22  —  (l/aii)a2iai2. 

Let  o-ij  ^  ^ij  ]*  Then 

a22  —  <^22  ~  (l/an)a2iai2. 


But,  by  assumption, 


detfA)  =  011022  —  021012  7^  0. 


Therefore, 


022  =  (l/aii)det(A)  7^  0. 


Appendix  4.c. 


(4.47) 


(4.48) 


(4.49) 


Proof  of  Theorem  4.5:  Select  any  arbitrary  vector  v  >  0.  Then,  since  A~^  >  0 
and  it  is  regular,  u  :=  A~*v  >  0.  Next,  let  A  €  [  A,  A  ],  obviously  A  <  A  <  A. 
Then, 

A"*A<I<A"^A.  (4.50) 

Therefore,  B  :=  A“^A  satisfies  B  <  I  and  Bu  =  A“^Au  >  A~^Au  =  A~^v  >  0. 
This  implies  that  B  is  an  Af-matrix  and  both  B  and  A  =  AB  are  regular. 

Now,  since  A“^  =  B~^A“*  >  0, 

A-^<A-^<A"V  (4.51) 

where  the  equality  holds  for  A  =  A  and  A  =  A.  B 
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Appendix  4.d. 


Proof  of  Theorem  4.8:  If  A  is  a  nonsingular  Af-matrix  then  there  exists  a  positive 
definite  W  such  that 

AW  +  WA”^ 


is  positive  definite.  Since  a  solution  W  exists,  then  it  is  clear  that  (I  0A^  +  A'^SI) 
is  a  nonsingular  matrix.  To  prove  the  result,  it  suffices  to  show  that  if  x  7^  0;  y  = 
Ax,  then  for  some  subscript  i,  r,  ^  0,  and  x,y,  >  0. 

For  simplicity  of  notation,  let  B  =  A^.  The  explicit  form  ofi5=(l0B  4-  B0l) 
is  given  by 


B  = 


BO  ¬ 

•  0 

bi2l 

•••  bynll 

OB 

•  0 

+ 

^21 1 

b22l 

•••  b2nl 

•  0 

•  0 

B. 

.bnll 

bn2l 

bnn^  - 

■  B  -h  61 1 1 

b\2 

I 

1 

621 1 

B  -f  622 1 

•  • 

b2nl 

. 

bn2 

I 

••  B-h5„„I. 

(4.52) 


By  assumption,  the  diagonal  elements  of  B  eire  positive,  therefore,  the  diagonal  ele¬ 
ments  of  B  eire  also  positive. 

2 

Now,  selecting  x  =  ^,,  the  i-th  column  of  the  identity  matrix  I  G  IR"  ,  it  is  easy 
to  see  that 


y  =  Bx  =  bi 


(4.53) 


where  bi  is  the  first  column  of  B.  Clearly,  for  i  =  1,  .ti  >  0  (=  1);  and  yi  =  2fcii 
(>  0  by  assumption).  Therefore  xij/i  >0,  proving  the  result.  I 
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5.  INTERVAL  ROUTH-HURWITZ  ARRAY 


AND  STABILIZATION  OF  UNCERTAIN  SYSTEMS 

5.1.  Introduction 

Analysis  and  design  of  systems  whose  parameters  are  known  to  lie  within  a  range 
(rather  than  having  an  exact  value)  has  received  considerable  attention  over  the  last 
decade.  Various  analysis  and  design  techniques  used  by  control  systems  engineers 
are  essentially  meant  for  application  to  a  “nominal”  model.  The  resulting  design  is 
said  to  be  robust  if  the  system  performs  within  acceptable  limits  in  the  face  of  sig¬ 
nificant  parameters  variations  emd  model  uncertainties.  There  is  a  huge  amount  of 
literature  that  addresses  the  problems  of  stability  and  stabilization  of  interval  poly¬ 
nomials.  The  first  paper  that  treated  the  problem  from  a  true  interval  point  of  view 
was  the  landmark  paper  by  Kharitonov  [22]  regarding  stability  of  polynomials  whose 
coefficients  are  independent  intervals.  Later  several  researchers  extended  the  result 
to  more  specialized  interval  polynomials,  where  the  coefficients  could  have  linearly 
coupled  uncertainties.  Notably,  the  works  of  Barmish  [4],  Keel  et  ai,  [21],  Saeki  [36] 
and  Zhou  and  Khargonekar  [45],  as  well  as  the  references  therein  provide  a  list  of 
various  efforts  in  this  direction. 

In  this  section,  we  define  and  develop  the  interval  Routh- Hurwitz  array  [25]  for  a 
transfer  function  whose  denominator  is  a  polynomial  with  interval  coefficients.  Based 
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on  the  interval  Routh-Hurwitz  {'RJi)  array,  we  will  present  a  procedure  for  finding  a 
state  feedback  vector  that  will  stabilize  the  entire  interv'al  system. 

5.2.  Interval  Routh-Hurwitz  Arrays 


In  this  Section,  we  will  develop  the  principles  on  which  the  interval  'R.'H  arrays  are 
based.  The  results  regarding  the  'R.'H  arrays  for  the  nominal  systems  are  extended  to 
the  case  of  polynomials  with  interval  coefficients. 

5.2.1.  PRELIMINARIES 


The  starting  point  for  the  subsequent  analysis  in  this  chapter  is  the  description  of  a 
linear  system  (single  input,  single  output)  is  in  a  differential  equation  of  the  following 
form: 


d" 


.  d”~^ 

[  d„_i,  a„_i  +  [  an-2,  d„_2  + 


im  —  l 


dV 


+  Cm-1 


d<"— 2 


+  •  •  •  +  C,n-2 


'dV 

d"* 


f  ^  ^  T  «  , 

+  t  ^  d? 


-2 


-I- 


+  C.JJ+CO 


u(f) 


y(<) 

(5.1) 


whose  transfer  function  is  given  by 

-f  H - hCiS4-Co 

G{s)  =  - p - ; - p - - - p - : — p - p - : — p-.  (0.2) 

s"  +  [  a„_i,  a„_i  ]s"-i  H - +  [  dp ,  op  Js  +  [  ao,  oq  j 

It  should  be  pointed  out  that  the  coefficient  of  the  highest  derivative  oi  y(t)  is  assumed 
to  be  [1,1],  i.e.,  the  denominator  polynomial  of  the  system  (5.2)  is  rnonic.  This 
assumption  is  not  necessary  for  developing  the  interval  IZH  array;  however,  it  is 
required  for  stabilization  of  (5.1)  by  state  feedback.  Note  that  such  an  assumption 
is  not  necessary  for  nominal  systems  because  the  polynomial  can  always  be  made 
monic  by  dividing  it  with  the  coefficient  of  the  highest  power  of  the  derivative.  This 
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operation  cannot  be  performed  on  the  interval  polynomials  because  a  multiplicative 
inverse  is  not  defined  for  interveil  arithmetic.  It  will  be  assumed  that  the  given  system 
is  controllable  for  all  variations  in  the  parameters. 

It  is  well  known  that  any  proper  scalar  rational  function  G(s)  =  Q(s)/P(s), 
where  P(s)  is  a  monic  polynomial,  can  be  written  in  a  controllable  but  not  necessarily 
observable  state  space  realization  given  by  [19]; 


■  0 

1 

0 

0  • 

‘O' 

0 

0 

0 

0 

0 

x(f)  = 

0 

0 

0 

1 

x(t)  -1- 

0 

.-Po 

-Pi  • 

1 

3 

1 

—Pn-1  . 

.1. 

y(t)  =  [90  92  •••  9n-2  9„-i  ]x(t)  +  du(t) 


(5.3) 


(5.4) 


where  pi  and  axe  the  coefficients  of  the  f-th  power  of  s  in  the  denominator  and  the 
numerator  polynomials,  respectively. 

Typically,  an  interval  7Z7i  array  for  a  scalar  system  (e.g.,  5-th  order  system)  with 
the  denominator  polynomial  having  intervrd  coefficients  as  shown  below 


d[s]  —  [1,  l]s®  -1-  [  02,1,  02,1  I'S'*  +  [  01,2,  di,2 

+  [  02,2,  02,2  1-S^  -I-  [  di,3,  01^3  -|-  [  02,3,  02,3  ],  (5.5) 


will  have  the  following  form: 


s'* 

[  02,1,  02,1  ] 

[  di,2,  di,2  ]  [  di,3,  di,3  ] 

[  02,3,  02,3  I  [  02,3,  02,3  ] 

s^ 

[  03,1 ,  03,1  I 

[  03,2,  03,2  ] 

[  04,1 ,  04,1  ] 

[  04,2,  04,2  I 

s“ 

[  05,1 ,  05,1  ] 

1  06,1  ,  06,1  ] 

where  d,  j  and  d,j  denote  the  lower  and  upper  bounds  of  various  interval  coefficients 
obtained  following  the  usual  rules  of  setting  up  tlu'  PTi  array  for  polynomials  with 
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constant  coefficients.  Specifically. 


[  ®3,1, 

«3,1  ]  = 

=  [[02,1, 

02,1 

][  Ol,2, 

di,2 

02,2, 

02,2  1]  /[  02,1  ,  02,1  1 

[  ^3,2, 

03,2  ]  = 

=  [[«2,1, 

02,1 

][  Ol,3, 

di  ,3 

]-[l,l][ 

02,3- 

02,3  ]]  /[  02,1 ,  02,1  ] 

[  ^4,1, 

04,1  ]  = 

=  [[«3,1, 

03,1 

][  02,2, 

02,2 

1  —  [  02,1 , 

02,1 

][  03,2,  03,2  ]]  /[  03,1  ,  03,1  ] 

[  «4,2, 

0> 

II 

=  [[«3,1, 

03,1 

][  ®2,3, 

02,3 

]]/[  5.3,1, 

03,1 

] 

[  Qs.i, 

05,1  ]  = 

=  [[04,1, 

04,1 

][  03,2, 

03,2 

]  “  [  53,1, 

03,1 

][  04,2,  04,2  ]]  /[  04,1 ,  04,1  ] 

[  °6,1, 

06,1  1  = 

=  [[«5,1, 

05,1 

][  ^4,2, 

04,2 

]]/[  55,1, 

05,1 

]  (5.6) 

The  only  modification  in  the  above  array  compared  to  the  'R.'H  arrays  for  nominal 
polynomials  is  that  to  compute  the  elements  of  various  rows,  we  use  the  principles  of 
interval  arithmetic  following  the  operations  reviewed  in  CHAPTER  2. 

It  was  shown  by  Krishnamurthy  [23]  that  for  a  stable  system,  the  elements  of 
any  row  of  the  'RfH  array  are  of  the  same  sign.  Further,  for  the  given  system,  if  any 
element  of  a  particular  row  is  found  to  be  of  different  sign,  it  can  be  concluded  that 
the  system  is  unstable  and  has  one  or  more  poles  in  the  right  half  plane.  In  the  next 
few  paragraphs,  we  establish  a  similar  result  for  the  interval  system  by  examining  the 
elements  of  the  interval  TZH,  array. 

5.2.2.  IMPLICATIONS  OF  INTERVAL  KH  ARRAYS 

It  was  shown  in  [23],  that  given  a  stable  nominal  polynomial 

F{s)  =  +  ^21^"  ^  +  ai2.s”~^  -!-  •  ••  a,  j  >  0  (5.7) 

such  that  all  elements  of  the  first  columns  of  the  TZTi  array  are  positive,  it  can  be 
concluded  that 
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Theorem  5.1:  For  a  stable  polynomial,  the  necessary  and  sufficient  condition  of  the 
'R.'H  criterion  implies  that  all  the  elements  of  the  TZ.'H  array  are  of  the  same  sign. 

Proof:  See  Appendix  5. A. 

For  an  n-th  order  polynomial,  the  following  'R'H  array  can  be  set  up: 


Ol.l  01,2  Ol,3  •  •  • 

02,1  02,2  02,3  ■  •  •  0.2, im 

^n-2 

03,1  03,2  03,3  •  •  •  03,j„ 

Oi,l  Oi,2  0,,3  •  •  •  0{,j^ 

On-2,1  0„_2,2 

On-1,1  0„_i,2 

On,l 

On+1,1 

where  jm  =  [(^  +  3  —  *)/2J  is  the  largest  integer  value  of  the  column  with  a  non-zero 
entry  in  the  i-th  row.  Note  that  the  above  notation  is  for  both  odd  and  even  degree 
polynomials.  However,  the  results  outlined  below  will  be  for  odd  degree  polynomials 
only.  Parallel  results  can  be  easily  obtained  for  even  degree  polynomials. 

For  a  stable  polynomial,  the  following  three  properties  hold: 

Property  1.  All  elements  in  the  first  column  of  the  array  are  positive 

>  0,  i  =  1,2,  ...,n-fl  (5.8) 

Property  2.  The  number  of  coefficients  in  alternate  rows  is  reduced  by  one, 

aij  =0,  i  =  1,2,  ...,n-f- 1,  j  >  [{n Z  -  i)/2\  (5.9) 

where  [(■)]  denotes  the  integer  part  of  the  number. 
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Property  3.  The  constant  term  appears  as  the  last  element  of  alternate  rows,  i.  c 


an+3-2k,k  =  constant  term,  k  =  1, 2, ... ,[(?!  + 2)/2)J .  (5.10) 

To  extend  the  result  to  the  interval  'R.'H  arrays,  consider  the  following  odd  piow- 
ered  interval  monic  polynomial: 

F{s)  =  [1,  l]s"  +  [  021,  021  ]•«"“*+[  0,12,  dl2  d,-,j  ]  >  0  (5.11) 

The  corresponding  interval  TZTi  array  is  given  by 


3”  1 

[  02,1.  02,1  ] 

[  Ol,2,  ai,2  ] 

[  02,2,  02,2  ] 

[  Ol,3,  Oi,3  ]  • 

[  02,3,  02,3  ]  • 

'  ■  [  Ol,j,„  ,  Oi.j,„  ] 

••  [  02,i,„,  a2,_,„,  ] 

[  03,1,  03,1  ] 

[  03,2,  03,2  ] 

[  03,3,  03,3  ]  • 

[o3,>„,  03,j„,  ] 

^n-i+1 

[  Ot,i ,  Sj  j  ] 

[  01,2,  dl,2  ] 

[  0,-,3,  01,3  ]  • 

■  ■  [  Ol,_7'„,  ,  ^i,jm  ] 

[  On-2,1,  dn-2,1  ] 

[  On-2.2,  dn-2,2 

] 

• 

[  On-1,1,  ] 

[  On-1,2,  On-1,2 

] 

s' 

[  On,l ,  dn,l  ] 

[  On+1,1,  On  +  1,1  ] 

where  all  elements  of  the  zu'ray  have  been  calculated  using  the  principles  of  interval 
arithmetic. 

Unlike  the  TIH  array  for  nominal  polynomials,  the  PROPERTY  3  does  not  hold 
for  interval  HH  arrays.  However,  the  following  is  true 

[  On+1,1-  dn.(-i  j  ]  D  [  On  — 1,2,  d„_i^2  ]  5  [  On-3,.3,  ^n-3,3  ]  ^  5  [  ^2,jm  ,  ]• 

(5.12) 

This  can  be  shown  trivially  following  the  sequence  of  operations  performed  to  ob¬ 
tain  the  interval  TZH  array.  This  property  will  be  referred  to  as  PROPERTY  3. 
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Clearly  by  PROPERTY  1,  [a„4.i,i,  ]  is  a.  positive  interval,  therefore,  the  in 

tervals  [a„_i,2,  a„_i,2  ]  [0^-3, 3,  an-3,3  ]  ■••  02, ]  are  also  positive. 

Following  the  same  reasoning  as  we  did  for  the  case  of  nominal  polynomials,  we 
next  develop  the  implications  of  the  interval  'R.'H  array.  Obviously,  from  PROPERTY 
1,  the  {n  +  l)-th  and  n-th  rows  have  positive  intervals.  Also,  from  PROPERTIES  1,  3, 
the  interval  in  (n  —  l)-th  row  are  positive.  In  (n  —  2)-th  row,  to  show  that  the  second 
interval  is  positive,  we  proceed  as  follows.  By  definition, 

r.  ..  i  [  On-l,!,  On-l,!  ][  a„_2,2,  ^n-2,2  ]  —  [  0„_2,1,  dn-2,1  ][  dn-1,2  ] 

lan,i,a„,ij  — - r-c - : - - 

(5.13) 

performing  the  inverse  interval  operations,  we  have 

[  Qn,l»  Qn,l  ][  Qn-l,!,  dn-1,1  ]  +  [  an-2,l>  Qn-2.1  ][  Qn-1,2>  Qn-1,2  ] 

[  Qn-1,1,  On-l,!  ] 

2  [  dn-2,2i  dn-2,2  ]•  (5.14) 

Each  interval  of  the  expression  on  the  left-hand  side  of  (5.14)  is  positive,  therefore. 
[  o„_2,2,  On-2,2  ]  is  Contained  in  a  positive  interval. 

In  the  (u  —  3)-th  row,  there  are  three  elements,  from  PROPERTIES  1,3,  the  first 
and  the  last  elements  are  positive  intervals.  To  show  that  the  second  element  is  also 

positive,  we  know  that 
[  ®n  — 1 ,1 1  dn  — 1 ,1  ] 

_  [  dn-2,1  ,  dn-2,1  ][  dn-3,2,  dn_3_2  ]  —  [  dn-3,1 ;  Qn-3,1  ][  dn-2,2)  Qn-2,2  ] 

[  dn-2,1,  d„_2,l  ] 

(5.15) 

From  where,  performing  the  inverse  operations,  we  get 

[  Qn-l,n  dn_i^l  ][  an_2j  ,  dn-2,1  ]  4-  [  Qn-3,1 1  ][  dn_2,2)  Qn-2,2  ] 

d„_2,I ,  d„_2,i  ] 
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(5.16) 


5  [  «n-3,2<  an-.3,2  ] 

Again,  since  each  element  of  the  on  the  right  hand  side  of  (5.16)  is  positive,  an-i,2 
also  positive. 

Continuing  in  this  manner,  we  show  that  the  elements  of  the  z-th  row  are  posirive 
intervals.  The  first  element  is  positive  by  PROPERTY  1.  The  7-th  element  of  the  ?-th 
row  is  given  by 

r  -  -  1  [  di+i,!,  di+i^t  [  d,7i ,  ] 

[ai+2,j-i,  ai+2,j-i  1  =  - j— - i - •  (5-1  0 

I  «i+l,0  ^1+1,1  J 

To  get  aji  inclusion  for  [  d,  j,  aij  ]  from  (5.17),  we  perform  the  inverse  operation  to 


[Q|-H,1i  qt+1,1  ][Q«+2,j-li  ai+2,j-l  ]  +  Qj.l  ][Qi+l.;7  Q«+l,j] 

[  Oj+i,!,  d<+i,i 


2  [  di.j ,  Qt.j  ]• 

(5.18) 


Since  all  of  the  intervals  on  the  left-hand  side  are  either  positive  by  PROPERTY  1  or 
have  been  proved  to  be  positive  in  the  previous  cycle,  the  interval  on  the  right  hand 
side  in  (5.18)  is  positive. 

Further,  for  even  z,  the  last  element  ]  is  included  in  the  interval 

[  ]  and  therefore  it  is  positive.  If  ?'  is  odd,  then 


^i  +  2,jm^  ^i  +  2,jrr,  Q«-H,l  ]  +  [  Qlfl.jVn,  ][  0»,1  ] 


Oi-l-i.o  ai-i-i.i 


5  ^  J- 

(5.19) 


which  is  a  positive  interval. 


Therefore  for  a  stable  interval  polynomial  with  positive  first  column  of  the  interval 
TZH  array,  every  element  of  the  array  is  positive.  Similar  results  can  be  obtained  when 
n  is  an  even  ordered  interval  polynomial. 


r,o 


It  should  be  pointed  out  that  unlike  the  classical  'RfH  array,  the  interval  'R'H 
array  provide  only  a  sufficient  condition  for  stability.  This  is  due  to  the  fact  that  in 
performing  the  interval  operations,  the  intervals  are  implicitly  overboimded. 

Example  5.1;  To  illustrate  the  above  point,  consider  the  polynomial  given  by 
F{s)  =  S®  +  [  021,  021  ]s^  +  [  012,  012  W  +  1- 

where  [  di2,  di2  ]  =  [4,6]  and  [  021,  021  ]  =  [2,3].  Now,  treating  the  two  intervals  as 

variable  coefficients  r2i  and  ri2,  we  get  the  following  nominal  'JtH  array: 

1  ri2 

r2i  1 


On  substituting  the  values  of  ri2  and  r2i  for  maximum  allowable  variation  in  the  two 
uncertain  coefficients,  we  get  the  following  array; 


1 

[4,6] 

[2,3] 

1 

[3.5,5.6667] 

0 

5° 

1 

0 

On  the  contrary,  if  the  intervals  were  allowed  to  “grow”  &s  the  array  was  computed, 
the  corresponding  interval  TZTi  array  will  be 


1 

14,6] 

12,3) 

1 

[2.3333,8.5] 

0 

[0.2745, 3.6429] 

0 

It  is  clear  from  the  above  example,  that  using  interval  arithmetic,  the  results  will  be 
over  bounded. 
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5.2.3.  REMOVAL  OF  CONSERVATISM 


The  example  (EXAMPLE  5.1)  presented  in  the  previous  section  provides  a  possible 
way  to  remove  the  conservatism  in  the  analysis  procedure.  Next,  we  formalize  one 
possible  way  for  removing  the  conservatism  in  developing  the  interval  'RfH  arrays. 

To  get  an  intuitive  understanding  of  the  procedure,  let  us  consider  the  element 
X31  in  the  first  column  of  .s^  row  from  EXAMPLE  5.1  given  by 


2^31  = 


ri2r2i  -  1 

^21 


(5.20) 


Notice  that  using  interval  arithmetic,  113  is  determined  as 

ri2r2i  -  1 


^13  ]  — 


r2i 

[4,6][2,3]  -  [1,1] 


12,31 

=  17,17111/3,1/2] 


(2.3333,8.5]. 


The  element  r2i  occurs  both  in  the  numerator  and  the  denominator.  The  interval 
inclusion  obtained  above  does  not  account  for  this  fact.  It  treats  each  interval  as 
independent  of  the  other  intervals  and  in  the  process  overbounds  the  interval. 

If,  on  the  other  hand,  we  account  for  the  fact  that  r2i  can  only  assume  one 
value  both  in  the  numerator  as  well  as  the  denominator,  then  the  interval  [  £13,  f  13  ] 
(=  [3.5, 5.6667])  so  obtained  is  the  maximum  variation  in  light  of  multiple  occurrence 
of  r2i .  This  fact  can  be  utilized  to  remove  the  conservatism  from  the  interval  'R'H 
eurrays. 

To  enable  us  to  obtain  such  minimixm  and  maximum  values  for  each  interval  in  tlu' 
first  column  of  the  interval  IZH  array,  we  need  a  numerical  algorithm  that  computes 
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the  extrema.  Fortunately,  Moore  [27]  and  Skelboe  [39]  provide  the  necessary  details 
for  developing  an  efficient  algorithm  for  obtaining  the  extrema.  A  detailed  account 
of  these  techniques  is  presented  by  Ratschek  and  Rokne  in  [34].  In  the  sequel,  we 
present  the  highlights  of  the  algorithm.  The  notation  that  we  have  adopted  is  similar 
to  the  one  in  [34]. 

The  problem  that  we  need  to  address  is  that  of  the  global  unconstrained  opti¬ 
mization  which  can  be  stated  as  follows:  Let  IR  be  the  set  of  reals,  let  X  C  IR"  be 
a  compact  right  parallelepiped  parallel  to  the  axes  (denoted  by  a  J),  /  :  X  — ♦  IR  any 
function  and  |/(X)  the  range  of  the  function  /  over  X,  i.e. 

B/(X)  =  {fix)  :  X  G  X}.  (5.21) 

The  global  minimum  (if  it  exists  in  IR),  defined  as  infimum  (|]/(X))  is  denoted  by  /. 
The  global  maximum  (if  it  exists  in  IR),  defined  as  supremum  ([[/(X))  is  denoted  by  /. 
Note  that  infimum  and  supremum  are  being  computed  instead  of  minimum  and  max¬ 
imum  because  of  lack  of  assumption  of  continuity  on  the  function.  Mathematically, 
the  problem  defined  above  can  be  stated  as 

minimize  /(x)  subject  to  x  G  X  (5.22a) 

maximize  /(x)  subject  to  x  G  X.  (5.226) 

In  the  following  paragraphs,  only  minimization  will  be  considered.  The  procedure  for 
maximization  is  identical  except  for  a  change  in  sign  of  /. 

Given  the  bounded  function  /,  the  domain  [X,  X  ]  G  HIR"  and  an  initial  choice  of 
inclusion  fianction  F  :  n]R(X)  — ►  HIR  of  /(x),  the  algorithm  works  by  splitting  up  the 
domain  X  in  each  step,  into  subboxes  of  not  necessarily  the  same  size.  The  search  for 
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/  is  done  in  the  subboxes.  A  direct  search  would  be  computationally  very  demanding, 
hence  a  branching  principle  is  used.  At  each  iteration,  the  search  is  continued  in  a 
box  Y,  where  /  has  the  smallest  lower  bound  y,  because  the  chances  of  finding  /  are 
the  best  in  this  box. 

Algorithm  5.1:  GLOBAL  UNCONSTRAINED  OPTIMIZATION 

Set  Y  :=  X 
Calculate  i^(Y) 

Set  y  :=  min  F{Y) 

Initialize  list  C  :=  ((Y,y)) 

Iterate: 

Select  a  coordinate  direction  k  parallel  to  which  Y  =  Yj  x  •  •  ■  x  Y„  has  an  edge 
of  maximum  length  i.e.  fc  €  {i  :  w{Y)  =  u)(Yj)} 

Bisect  Y  normal  to  direction  k  obtaining  boxes  Vj  and  V2  such  that  Y  =  Vj 
U  V2 

Calculate  F(Vi)  and  F(Vi) 

Set  Vi  =  Fib{Vi),  i  =  l,2 
Delete  (Y,y)  from  the  list  C 

Include  (Vi,V])  and  (V2,V2)  in  the  list  such  that  second  members  of  all 
pairs  of  list  do  not  decrease 
Denote  the  first  peiir  of  the  list  by  (Y,y) 

If  w(F{Y k})  <  e  then  terminate  else  iterate 

The  detailed  implementation  of  the  algorithm  as  well  as  further  details  regarding 
the  convergence  and  termination  criteria  can  be  found  in  [34].  It  suffices  to  say  that  the 
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necessaxy  tools  for  obtaining  the  correct  bounds  on  the  intervals  in  the  first  column  of 
the  interval  TtTi  array  are  available.  It  is  also  worth  mentioning  that  several  software 
packages  that  implement  machine  interval  arithmetic  exist.  Some  of  the  commonly 
used  ones  are  TRIPLEX-ALGOL-60,  PASCAL-SC,  FORTRAN-SC,  ACRITH  and 
ARITHMOS. 

Using  above  algorithm,  we  can  make  the  condition  on  the  interval  Ttl-t  array  both 
necessary  and  sufficient.  The  modified  procedure  would  require  the  determination  of 
the  elements  of  the  first  column  of  the  array  symbolically  and  then  application  of  above 
algorithm  will  provide  us  with  the  correct  lower  and  upper  bounds  of  each  element 
in  the  first  column  of  the  interval  'R/H  array.  For  symbolic  computation  software 
packages  like  MathemaTICA  can  be  used. 

5.3.  Stability  of  Systems  with  Dependent  Uncertainties 

Unlike  systems  described  in  earlier  section,  where  each  coefficient  has  variations  that 
are  totally  independent  of  the  rest,  many  systems  exhibit  dependent  uncertainties  [4], 
[21],  [36]  [45].  Specifically,  each  coefficient  may  have  several  uncertain  parameters. 
Further,  these  peirameters  may  exhibit  some  functional  behavior.  In  this  section,  we 
will  extend  the  results  of  SECTION  5.2.3  to  systems  with  dependent  uncertainties.  It 
will  be  shown  that  in  view  of  the  previous  section,  the  stability  of  interval  polynomi¬ 
als  with  dependent  uncertainties  can  be  determined  in  a  simple  manner.  Note  that 
several  papers  deal  with  the  above  issue,  however,  the  techniques  tend  to  be  either  too 
complicated  or  computationally  extremely  demanding.  In  the  next  few  paragraphs 
we  show  that  the  procedure  outlined  in  the  previous  section  could  be  easily  applied 
to  study  the  stability  properties  of  polynomials  whose  coefficients  exhibit  linear  and 
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non-lineax  dependent  uncertainties.  Formal  proofs  of  the  material  presented  in  the 
sequel  have  not  been  developed.  Hence,  only  examples  have  been  used  to  illustrate 
the  concepts. 

Let  us  denote  the  characteristic  polynomial  of  the  systems  with  coupled  uncer¬ 
tainties  as 

n  — 1 

d[s]  =  s"  +  v;[/.(r)]^-  (5.23) 

«=o 

where  /t(r)  denotes  the  coefficient  of  the  i-th  power  of  s  and  r  is  the  vector  of  uncertain 
terms. 

We  discuss  three  possible  ways  to  treat  the  situation  discussed  above. 

1.  The  first  approach  could  be  to  overbound  the  various  coefficients  using  interval 
arithmetic  directly,  and  use  interval  7l7i  criterion.  As  expected  this  would  lead  to 
conservative  results. 

Example  5.2;  To  illustrate  this  point,  consider  the  polynomial  described  by 

d[s]  =  +  [ri+r2]s^  +  [rjr2]s  +  [fj], 

where,  rj  G  [0.74,0.76]  and  r2  €  [2.00,2.50].  The  polynomial  has  two  uncertain 
elements  which  appear  as  coupled  uncertainties  in  various  coefficients.  On  expanding 
the  polynomial,  the  numerical  values  in  various  coefficients  are  given  by 

d[s]  =  +  [2.74,3.26]^^  -I-  [1.48,1.9].s  +  [4,6.25]. 

Now,  using  interval  TZH  array,  we  have: 


1 

[1.48,1.90] 

[2.74,3.26] 

[4.00,6.25] 

[-0.8010,  .8007] 

— 

undefined 

0 

.%G 


Note  that  in  the  above  array,  the  first  element  in  the  third  row  includes  zero,  since 
the  elements  of  the  subsequent  rows  would  require  division  by  this  element,  we  cannot 
proceed  any  further.  Clearly  the  criterion  suggests  a  sign  change  in  the  first  column 
and  hence  instability. 

2.  The  second  possibility  is  to  verify  whether  the  expanded  polynomial  meets 
Kharitonov’s  criterion.  This  entails  checking  the  stability  of  the  following  four  poly¬ 
nomials: 

-I-  3.26s^  +  1.48s  -I-  4.00 
+  2.743^  -1-  1.90s  -h  6.25 
s^  -f-  3.26s2  -I-  1.90s  -f-  4.00 
s^  -I-  2.74s2  -t-  1.48s  -t-  6.25. 

While  the  first  eind  third  polynomials  are  stable,  the  second  and  the  fourth  polynomials 
are  unstable,  hence  use  of  Kharitonov  criterion  on  the  expanded  polynomied  also 
indicates  instability. 

3.  Finally,  the  third  approach  would  be  to  use  the  Global  Unconstrained  Optimization 
on  the  elements  of  the  symbolic  TT/H  array  and  infer  stability  from  there.  The  array 

takes  the  following  form: 

s^  1  rir2 

s^  ri  -t-  r2 _  r| 

^2  (^1  -  7^)  ^ 

s«  1 

Now,  using  unconstrained  optimization,  the  first  element  in  the  third  row  is  found 
to  be  [0.0201,  0.0884],  which  is  positive  and  hence  the  polynomial  with  dependent 
coefficients  is  stable. 
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The  technique  is  simple  to  apply  because  the  optimization  technique  is  commer¬ 
cially  available  in  packages  that  implement  interv'al  arithmetic. 

5.4.  Feedback  Stabilization  of  Uncertain  Systems 

This  section  presents  results  on  the  use  of  the  techniques  developed  above  in 
robust  synthesis  of  feedback  control  systems.  The  results  presented  in  this  section  are 
for  single  input,  single  output  systems  only. 

5.4.1.  SYSTEMS  WITH  INDEPENDENT  UNCERTAINTIES 

In  this  section,  we  will  show  that  by  using  the  interval  'R.'H  arrays,  we  can  determine 
a  state  feedback  vector  k  (when  it  exists)  such  that  the  closed  loop  “interval’’  system 
is  always  stable.  For  the  purpose  of  illustration,  we  will  develop  all  results  for  a  fifth 
order  system.  The  final  results  will  be  stated  for  the  general  order.  It  is  assumed  that 
the  uncertainties  in  various  coefficients  are  mutually  independent.  Further  assume 
that  the  system  is  in  the  following  canonical  form: 

dx(t)  [0  I  In-i]  ,,,  ,  ^  ,,, 

+  *^“(0  (5.24) 

yit)  =  cx(0  (5.25) 

where 

[^]^  =[  ~[  02, 3i  5:2,3  ]  —  (  Oi,3,  ai,3  ]  —  [  02,2 >  02,2  ]  —  [  Oi,2,  01,2  ]  —[  02,1  .  02,1  ]  ] 

=  [  0  0  0  0  1  ] 
c  =  [  C2,3  Cl  3  C2.2  <-1,2  C2,|  ] 

Define  a  state  feedback  vector  =  [1-2,3  A'1,3  ■  •  •  ^2,1].  Then  the  state  matrix 
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of  the  closed  loop  state  space  system  under  the  feedback  law 


u{t)  =  v{t)  -  k^x(<) 


(5.26) 


becomes 


fO  I  I„-il 
.  [a(k)]^  . 


where  [a(k)]  is  given  by 


[a(k)]  = 


■  — [  02,3,  02,3  ]  —  ^2,3  ' 

—  [  ai,3,  di,3  ]  —  fci,3 

—  [  02,2,  02,2  ]  —  ^2,2 
—  [  01,2,  di,2  ]  —  ^1,2 

.— [  02,1,  02,1  ]  —  ^2,1  - 


(5.27) 


Since  Ad  is  in  the  companion  form,  its  characteristic  polynomial  is  given  by 


^c[^]  =  [1,  Ij-S*  +  [o2,l  +  ^2, 1,02,1  +  A:2, IIS'* 

+  [di,2  +  ^1,2,01,2  +  A:i,2]s^  +  [02,2  +  ^2, 2, 02,2  +  ^2,2]'S^ 

+  [dl,3  +  ^1,3)01,3  +  fcl,3]'S*  +  [02,3  +  ^2, 3, 02,3  +  ^2, si's”,  (5.28) 


and  the  corresponding  interval  TZH  array  is  modified  to 


s 

s 

s 

s 

s 

s 


5 

4 

J 

2 

1 

0 


[1,1] 

[d2,i(k),a2,i(k)] 

[d3,i(k),a3,i(k)j 

[04,1  (k),  04,1  (k)] 
[05,1  (k),  05,1  (k)] 

[d6,i(k),a6,i(k)] 


[di,2(k),ai,2(k)] 

[q2,2(k),d2,2(^)] 

[d3,2(k),d3,2(k)] 

[d4,2(k),a4,2(k)] 


[di,3(k),ai,3(k)] 

[d2,3(k),a2,3(k)] 


where  a,,j(k)  and  d,_j(k)  are  now  functions  of  the  elements  of  the  state  feedback 
vector  k.  Now  extending  the  property  of  classical  7Z7i  arrays  to  interval  7Z7{  arrays, 
we  have  the  following  result: 
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Theorem  5.2  A  sufficient  condition  for  all  roots  of  the  interval  polynomial  to  lie 
in  the  left  half  plane  is  that  the  intervals  [  a,__,(k),  o;  j  (k)  ],  ?  =  1, 2, . . . ,  n  -H  1.  V  j , 
of  the  interval  'RfH  array  have  the  same  sign. 

To  solve  for  the  elements  of  k,  assuming  that  the  coefficient  of  the  highest  power 
of  s  is  [1,1],  we  set  up  the  following  n  non-linear  algebraic  inequalities  in  n  unknowns; 

fl2,i(k)  >  0 

03,1  (k)  >  0 

a„+,,i(k)  >  0.  (5.29) 

It  is  fairly  straightforward  to  see  that  if  a  solution  to  above  simultaneous  inequalities 
exists,  it  will  provide  us  with  an  n  dimensional  coefficient  space,  any  element  from 
which  would  ensure  the  stability  of  the  interval  system. 

Note  that  by  the  definition  of  an  interval,  it  is  assured  that  the  first  element  of 
the  interval  is  smaller  than  the  second.  Hence,  if  the  lower  extreme  of  the  intervals  in 
the  first  columns  have  positive  sign,  then  the  upper  extremes  are  guaremteed  to  have 
the  same  sign.  Further,  from  the  results  in  SECTION  3.2,  every  element  of  each  row 
will  have  the  same  sign,  satisfying  the  condition  in  [23].  The  above  design  procedure 
is  illustrated  by  means  of  a  3rd  order  system. 

Example  5.3:  Assume  that  the  denominator  polynomial  of  the  interval  system  is 
4s]  =  [l,l].s'  -f-  [-2,]]s"  +  [2,3]s  +  [0,2]. 

Clearly,  the  interval  polynomial  is  not  stable.  We  need  to  find  a  state  feedback  vector 
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=  [ko  ki  fca]  such  that 


+  [— 2  +  ^2, 1  +  ^’2]s^ 

+  [2  +  ^1,3  +  fci].s  +  [0  4-fco?2  +  fco] 

is  stable  for  the  entire  interval.  Using  the  interval  IZTi  array,  and  a  fair  amount 
of  painful  symbolic  manipulation,  one  gets  the  following  set  of  non-linear  algebraic 
inequalities: 


— 2  -|-  k2  >0 
(2 -|- fci)(— 2 ^,’2)  —  (2 -|- fco)  >0 
A:o(-2  + A:2)((2-|-iti)(-2-l-fc2)  ~  (2  +  ko))  >0 

that  define  the  region  from  which  the  parameters  of  k  can  be  obtained  to  stabilize 
the  system.  One  possible  stabilizing  feedback  vector  is  obtained  as  k  =  [1  2  3].  The 
corresponding  closed  loop  characteristic  polynomial  is  given  by 

dc[s]  =  [l,l]s*  -|-  [l,4]s^  [4,5]s  +  [1>3]- 

It  can  be  readily  verified  that  dc[s]  is  stable  for  entire  range  of  interval  coefficients. 

Since  the  stability  condition  presented  above  is  only  sufficient,  clearly  the  design 
so  obtained  will  be  conservative  in  nature.  However,  the  usefulness  of  employing 
interval  arithmetic  is  self-evident  in  that  it  provides  us  with  a  solution  (when  it  exists) 
in  a  fairly  straightforward  manner. 

It  is  easy  to  see  that  the  above  results  can  bo  used  for  simultaneous  stabilization 
by  means  of  state  feedback.  Assume  that  we  arc  given  N  plants  with  characteristic 
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polynomials 


d^^](s)  =  .s"  + 

j=0 

(^[^1(5)  z=  5”  -|- 

i=0 


d^'^{s)  =  s”  +  (5.30) 

i=0 

Assume  further  that  the  order  of  each  plant  is  n.  Then,  we  can  easily  obtain  an 
interval  plant  defined  as 

n— 1 

d[s]  =  s"  -f-  (5.31) 

1=0 

where  =  min(a[^^ )  and  =  max(a['’^ ),  j  =  1, 2, . . . ,  N.  If  a  controller  that  can 
stabilize  the  above  interval  plant  can  be  found  then,  clearly,  it  will  stabilize  all  the 
plants  in  the  above  family. 

Appendix  s.a. 

Proof  of  Theorem  5.1:  Although  the  proof  appears  in  [23],  an  abbreviated  version 
of  the  proof  is  included  for  the  sake  of  completeness.  For  an  n-th  order  polynomial, 
the  following  IZH.  array  can  be  set  up: 


s" 

01,1 

02,1 

Ol,2 

02,2 

Ol,3  •  • 
02,3  •  •  ' 

'  <^l,>m 
02,7m 

O3,] 

03,2 

03,3  •  •  ' 

^n-.+  l 

Oi,l 

Oi,2 

Oi,3  •  •  ■ 

.s’* 

an-2,1 

0 

T  •  ■ 

.s« 

On-l,! 

On,l 

On+1,1 

On-1,2 
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where  jm  =  [(”  +  3  —  0/2j  is  the  largest  integer  value  of  the  column  with  a  non- zero 
entry  in  the  i-th  row. 

For  a  stable  polynomial,  the  following  three  properties  (repeated  here  for  the  sake 
of  completeness)  hold: 

Property  1.  All  elements  in  the  first  column  of  the  array  are  positive 

Oi.i  >  0,  i  =  1,2,  ...,n -f- 1  (5.32) 

Property  2.  The  number  of  coefficients  in  2ilternate  rows  is  reduced  by  one, 

aij  =  0,  i  =  l,2,...,n-|-l,  j  >  [(n-h3-0/2J  (5.33) 

where  [(•)]  denotes  the  integer  part  of  the  number. 

Property  3.  The  constant  term  appears  as  the  last  element  of  alternate  rows,  i.e. 

an+3-2fc,fe  =  constant  term,  k  =  1, 2, . . . ,  [(n -I- 2)/2)J .  (5.34) 

Now,  from  PROPERTY  1,  the  (n  -|-  l)-th  and  n-th  rows  are  clearly  positive.  From 
Properties  1,  3,  the  elements  in  (n  —  l)'th  row  are  positive.  In  (n  —  2)-th  row,  to 
show  that  the  second  element  is  positive, 

an-2,2  =  (5.35) 

On-1,1 

Since  each  element  of  the  elements  on  the  right  hand  side  of  (5.35)  is  positive,  a„_2,2 
is  also  positive. 

In  the  (n-3)-th  row,  there  are  three  elements,  from  PROPERTIES  1,3  the  first  and 
the  last  elements  are  positive.  To  show  that  the  the  second  element  is  also  positive, 

On-:!, 2  =  «n-l,l  +  ^’’«n-2,2-  (5.36) 

flri-2,1 
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Again,  since  each  element  of  the  elements  on  the  right  hand  side  of  (5.36)  is  positive 
Qn-3,2  is  3Jso  positive. 

Proceeding  in  a  similar  manner,  consider  the  elements  of  the  i-th  row.  The  first 
element  is  positive  by  PROPERTY  1. 


ai,2  —  +  - ^1+1,2 


(5.31 


where  ai+1,2  would  have  been  proved  to  be  positive  from  previous  iteration,  hence 
0^,2  is  positive.  Talcing  the  general  term 


a. 


Oi,l 

- Oi+l 


(5.38) 


is  positive  by  the  same  reasoning  as  above. 

For  even  t,  the  last  element  is  the  same  as  the  constant  term  and  therefore 
it  is  positive.  If  i  is  odd,  then 


-  + 


<^1+1.1 


®i+l  .in 


(5.39) 


Therefore,  aij^  is  also  positive.  Similar  results  hold  when  n  is  even  ordered. 
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6.  RECOMMENDATIONS  FOR  FUTURE  WORK 


The  results  presented  in  this  report  are  preliminary  findings  on  the  use  of  interval 
arithmetic  in  analysis  and  design  of  control  systems.  Considerable  amount  of  work 
needs  to  be  done  to  make  it  viable  for  applications  in  practical  systems. 

Some  of  the  issues  that  were  not  addressed  in  this  report  and  would  make  inter¬ 
esting  topics  for  extension  of  this  work  are: 

1.  For  the  state  space  representation  the  analysis  was  restricted  to  the  class  of 
systems  where  the  state  matrix  is  an  Al-matrix.  A  more  general  framework  that 
addresses  a  broader  class  of  systems  is  certainly  desirable.  This  would  encompass 
developing  necessary  and  sufficient  conditions  for  stability  of  less  restrictive  class 
of  matrices.  Several  results  (see  list  of  references)  in  this  direction  have  been 
reported,  however,  they  tend  to  be  computationally  very  demanding.  Additional 
work  is  required  to  make  these  methods  computationally  efficient. 

2.  Possibilities  of  the  use  of  symbolic  arithmetic  should  be  explored  to  reduce  the 
conservatism  experienced  in  analysis  techniques.  Symbolic  arithmetic  together 
with  Global  Unconstrained  Optimization  presented  in  SECTION  5.2  could  provide 
the  necessary  machinery  to  circumvent  the  problem  of  conservatism. 

3.  The  topic  of  multi-parameter  as  well  a.s  non-linear  uncertainties  can  be  dealt  with 
following  the  approach  presented  in  Sectioii  5.3.  Formal  methodologies  as  well  as 
mathematical  justification  needs  to  be  investigated. 
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4.  The  results  presented  for  state  feedback  design  are  restricted  to  single  input  single 
output  systems.  An  extremely  important  research  problem  could  be  to  explore 
the  possibility  of  extending  them  to  multi-input  systems. 

5.  In  practical  control  systems,  state  feedback  may  not  always  be  feasible,  hence 
investigation  of  observer  based  feedback  as  well  as  output  feedback  problems  are 
important.  One  possible  approach  to 
the  generalization  of  methods  that  use 


address  output  feedback  control  may  be 
Sylvester  equations  for  feedback  control. 
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